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This article presents practical numerical recipes for simulating high-temperature and non- 
equilibrium quantum spin systems that are continuously measured and controlled. The no- 
tion of a "spin system" is broadly conceived, in order to encompass macroscopic test masses 
as the limiting case of large- j spins. The simulation technique has three stages: first the 
deliberate introduction of noise into the simulation, then the conversion of that noise into an 
informatically equivalent continuous measurement and control process, and finally, projec- 
tion of the trajectory onto a Kahlerian state-space manifold having reduced dimensionality 
and possessing a Kahler potential of multilinear (i.e., product-sum) functional form. These 
state-spaces can be regarded as ruled algebraic varieties upon which a projective quantum 
model order reduction (QMOR) is performed. The Riemannian sectional curvature of ruled 
Kahlerian varieties is analyzed, and proved to be non-positive upon all sections that contain 
a rule. It is further shown that the class of ruled Kahlerian state-spaces includes the Slater 
determinant wave-functions of quantum chemistry as a special case, and that these Slater 
determinant manifolds have a Fubini-Study metric that is Kahler-Einstein; hence they are 
solitons under Ricci flow. It is suggested that these negative sectional curvature proper- 
ties geometrically account for the fidelity, efficiency, and robustness of projective trajectory 
simulation on ruled Kahlerian state-spaces. Some implications of trajectory compression 
for geometric quantum mechanics are discussed. The resulting simulation formalism is 
used to construct a positive P-representation for the thermal density matrix and to derive 
a quantum limit for force noise and measurement noise in monitoring both macroscopic 
and microscopic test-masses; this quantum noise limit is shown to be consistent with well- 
established quantum noise limits for linear amplifiers and for monitoring linear dynamical 
systems. Single-spin detection by magnetic resonance force microscopy (MRFM) is then 
simulated, and the data statistics are shown to be those of a random telegraph signal with 
additive white noise, to all orders, in excellent agreement with experimental results. Then 
a larger-scale spin-dust model is simulated, having no spatial symmetry and no spatial 
ordering; the high-fidelity projection of numerically computed quantum trajectories onto 
low-dimensionality Kahler state-space manifolds is demonstrated. Finally, the high-fidelity 
reconstruction of quantum trajectories from sparse random projections is demonstrated, 
the onset of Donoho-Stodden breakdown at the Candcs-Tao sparsity limit is observed, and 
methods for quantum state optimization by Dantzig selection are given. 
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1. Introduction 

This article describes practical recipes for the simulation of large-scale open quantum spin 
systems. Our overall objective is to enable the reader to design and implement practical 
quantum simulations, guided by an appreciation of the geometric, informatic, and algebraic 
principles that govern simulation accuracy, robustness, and efficiency. 

1.1. How does the Stern- Gerlach effect really work? 

This article had its origin in a question that Dan Rugar asked of us about five years ago: 
"How does the Stern-Gerlach Effect really work?" The word ^^really" is noteworthy because 
hundreds of articles and books on the Stern-Gerlach effect have been written since the orig- 
inal experiments in 1921 |88| l89l [90] . . . including articles by the authors |178[ I182j and 
by Dan Rugar |179] himself. Yet we were unable to find, within this large literature, an 
answer that was satisfactory in the context in which the question was asked, that circum- 
stance being the (ultimately successful) endeavor by Rugar's IBM research group to detect 
the magnetic moment of a single electron spin by magnetic resonance force microscopy 
(MRFM) [T70] . 

1.1.1. Constraints upon the analysis Quantum theory has a reputation for mystery. But 
as Peter Shor has remarked, "Interpretations of quantum mechanics, unlike Gods, are not 
jealous, and thus it is safe to believe in more than one at the same time." In particular, it 
is well known — and we will review the literature in this article — that advances in quantum 
information theory have provided Shor's principle with rigorous foundations. 

We will build upon these informatic foundations in answering Dan Rugar's question 
in accord with the following constraints: our analysis will be orthodox in its respect for 
established principles of quantum physics. It will be operational in the sense that all its 
predictions are traceable to explicitly hardware-based measurement processes. The analysis 
will be scalable to accommodate large-dimension quantum systems (such as the spins in 
protein molecules that are the ultimate targets of MRFM microscopy). The analysis will be 
reductive in the sense that the analysis will yield simple design rules that are in reasonable 
quantitative accord with the predictions of more accurate — but more complicated — large- 
scale numerical simulations. The analysis will be synoptic in the sense that when we are 
required to choose between equivalent analysis formalisms, a rationale for these choices will 
be provided, and the consequences of alternative choices noted. And finally, the analysis 
will be extensible — at least in principle — to the analysis and simulation of general quantum 
systems (such as spintronic devices, nanomechanical devices, and biomolecules). 

There are of course strong practical motivations for seeking to analyze quantum systems 
by methods that are orthodox, operational, scalable, reductive, synoptic, and extensible: 
these same attributes are essential to practical methods for analyzing large-scale classical 
systems [201 [m]. 

1.2. The feasibility of generic large-scale quantum simulation 

We did not begin our investigations with the idea that the numerical simulation of large-scale 
quantum spin systems was feasible. Indeed, we were under the opposite impression, based 
upon the no-simulation arguments of Feynman [67] in the early 1980s. These arguments 
have been widely — and usually uncritically — repeated in textbooks |146[ sec. 4.7]. But 
Feynman's arguments do not formally apply to noisy systems, and in the course of our 
analysis, it became apparent that this provides a loophole for developing efficient simulation 
algorithms. Furthermore, it became apparent that the class of noisy systems encompasses 
as a special case the low-temperature and strongly correlated systems that are studied in 
quantum chemistry and condensed matter physics. In the concluding section of this article 
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(Section 4.6), we will develop the point of view that any quantum state that has been in 
contact with a thermal reservoir is an algorithmically compressible object. 

This loophole helped us understand why — from an empirical point of view — simulation 
capabilities in quantum chemistry and condensed matter physics have been improving ex- 
ponentially in recent decades j80| 11491 1151] . The analysis and simulation methods that we 
will present in this article broadly define a geometric and quantum informatic program for 
sustaining this progress. 

1.2.1. The geometry of reduced- order state-spaces This article's mathematical methods are 
novel mainly in their focus upon the geometry of reduced-order quantum state-spaces. We 
will show that the quantum state-spaces that are most useful for large-scale simulation 
purposes generically have an algebraic structure that can be geometrically interpreted as 
a network of geodesic curves (rules) having nonpositive Gaussian curvature for all sections 
that contain a rule (Section |2|. We will see that these curvature properties are essential to 
the efficiency and robustness of model order reduction. 

1.2.2. The central role of covert measurements A technique that is central to our simula- 
tion recipes is to simulate all noise processes (including thermal baths) as equivalent covert 
measurement and control processes (Section [3|. From a quantum informatic point of view, 
covert quantum measurement processes act to quench high-order quantum correlations that 
otherwise would be infeasibly costly to compute and store (Section [4]). Thus the presence 
of noise can allow quantum simulations to evade the no-simulation arguments of Feynman 

1.2.3. Background assumed by the presentation No reader will be expert in all of the 
disciplines that our analysis and simulation recipes embody, which are (chiefly) quantum 
mechanics in both its physical and informatic aspects, the engineering theory of model 
order reduction (MOR) and dynamical control, and the mathematical tools and theorems 
of algebraic and differential geometry. Indeed, the writing this article has made the authors 
acutely aware of their own considerable deficiencies in all of these areas. Recognizing this, we 
will describe all aspects of our recipes at a level that is intended to be broadly comprehensible 
to nonspecialists. 

1.2.4. Overview of the analysis and simulation recipes We begin by surveying our chain 
of reasoning in its entirety. Figures [l}|4] concisely summarize the simulation recipes and 
their geometric basis. In a nutshell, the recipes embody the orthodox quantum formalism 
of Fig. [T] as translated into the practical numerical algorithm of Fig. |2] which is based 
upon the algebraic structures of Fig. [3j whose functionality depends upon the fundamental 
geometric concepts of Fig. |4| 

The following overview summarizes those aspects of the simulation recipes that are 
fundamental, multidisciplinary, or novel, and it also seeks to describe the embedding of 
these recipes within the larger literature. 

1.3. Overview of the formal simulation algorithm 

Formally, our simulation algorithms will be of the general quantum information-theoretic 
form that is summarized in Fig. [T] Steps A. 1-2 of the algorithm are adopted, without 
essential change, from the axioms of Nielsen and Chaung |146] . and our discussion will 
assume a background knowledge of quantum information theory at the levels of Chapters 2 
and 8 of their text. 

Step A. 3 of the simulation algorithm — projection of the quantum trajectory onto a state- 
space of reduced dimensionality — will be familiar to system engineers as projective model 



order reduction (as we will review in Sec. 1.4). We will also establish that projective MOR 
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Formal Simulation of a Quantum Trajectory 

Notation: let H he a Hamiltonian and let measurement 
and noise operations be {{M^^^, M^^^}; = 1, . . . , m} . 



Optional: let /C be a Kahler submanifold. For j?/'^) G /C, 
let {Pic)n project onto the tangent space of /C at lipn) ■ 

Evolve the trajectory {\iJji),\iJj2), ■ ■ ■ jltpn), ■ ■ ■} as follows: 

while the simulation is running 

/* Step A.l: evolve the closed system dynamics */ 

let iV'n+i) = ex.p{-iHSt/h) jV^^) 

/* Step A. 2: evolve the open system dynamics */ 

for each of k = 1,2, ... ,m pairs of operations 



end while the simulation is running 



Figure 1: Formal simulation by quantum model order reduction. 

Steps A. 1-2 summarize the formal theory of the simulation of quantum systems (see, e.g., 
Nielsen and Chuang |146| chs. 2 and 8]). Step A. 3 is a model order reduction of the 
Hilbert states \jpn) by projection onto a reduced-dimension Kahler manifold fC (see e.g. 
Rewiehski [164!). Equivalently, Step A. 3 may be viewed as a variety of Dirac-Frenkel vari- 
ational projection (see, e.g., |133| 1161] ). 



let = {^n+l\M^J M^^^\llJn+l) 



/* note : p|V) + pf_) = 

{with probability 
{otherwise) 



1 */ 



end for each of m pairs of operatic 
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is formally identical to a method that is familiar to physicists and chemists as a variational 
order reduction of Dirac-Frenkel-McLachlan type |59| 176] I136j (see also the recent references 

laHiiiMlllll]). 

For purposes of exposition, we define quantum model order reduction (QMOR) to be 
simply classical MOR extended to the complex state-space of quantum simulations. 

1.3.1. The operational approach to quantum simulation Our simulation recipes will adopt 
a strictly operational approach to measurement and control, in the sense that we will require 
that the only information stream used for purposes of communication and control is the 
stream of binary stochastic outcomes of the measurement operations of Step A. 2 of Fig. [T} 
Although it is not mathematically necessary, we will associate these binary outcomes with 
the classical "clicks" of physical measurement apparatuses, and we will develop a calibrated 
physical model of these clicks that will guide both our physical intuition and our simulation 
design. 

1.3.2. The embrace of quantum orthodoxy Because the binary "clicks" of measurement 
outcomes are all that we seek to simulate, our analysis will regard the state trajectories {ipi, 
■02; • • • 5 'V'n) • • • } wholly inaccessible for all purposes associated with measurement and 
control, which is to say, as inaccessible for engineering purposes. We will analyze quantum 
state trajectories only with the goal of tuning the simulation algorithms to compress the 
trajectories onto low-dimension manifolds. In practice, this will mean that we mainly care 
about the geometric properties of quantum trajectories; this will be the organizing theme 
of our analysis. 

In the course of our analysis we will confirm — mainly to check our algebraic manipu- 
lations — that several of the traditional quantum measurement short-cuts that deal directly 
with wave-functions {e.g., uncertainty principles, wave function collapse, quantum Zeno 
effects) yield the same results as our "clicks-only" reductive formalism. But our simulations 
will not use these short-cuts, and in particular, we will never simulate quantum measurement 
processes in terms of von Neumann-style projection operators. 

The resulting simulation formalism is wholly operational, and can be informally de- 
scribed as "ultra-orthodox." The operational approach will require some extra mathe- 
matical work — mainly in the area of stochastic analysis — but it will also yield some novel 
mathematical results, including a closed-form positive P-representation [153" of the thermal 
density matrix. We will derive this P-representation by methods that provably simulate 
finite-temperature baths. Thus the gain in practical simulation power will be worth the 
effort of the extra mathematical analysis. 

1.3.3. The unitary invariance of quantum operations Our analysis will focus considerable 
attention upon the sole mathematical invariance of the simulation algorithm of Fig.[T| which 
is a unitary invariance associated in the choice of the quantum operations M in Step A. 2. 
Our main mathematical discussion of this invariance will be in Section 3.2.1[ our main 



discussion of its causal aspects will be in Section 3.3.6, our main review of the literature 



will be in Section 3.3. 7[ and it will be central to the discussion of all the simulations that 



we present in Section |4j 

We will see that the short answer to the question "What is this unitary invariance all 
about?" is that (1) it ensures that measured quantities respect physical causality, and (2) it 
allows quantum simulations to be tuned for improved efficiency and fidelity. 

In preparation, we caution readers that what we will call "quantum operations" are 
known by a great many other names too, including Kraus operators, decomposition opera- 
tors and operation elements. These operations are discussed in textbooks by Nielsen and 
Chuang |146j . Alicki and Lendi [6], Carmichael |40| 141], Percival [ 152] , Breuer and Petruc- 
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cione [251 , and Peres |155j . These texts build upon the earher work of the mathematicians 
Stinespring [186J and Choi [H] and the physicists Kraus |122| 1123] . Davies [57 , and Lind- 
blad |131j . Shorter, reasonably self-contained discussions of open quantum systems can be 
found in articles by Peres and Terna |155] . Adler p], Rigo and Gisin |166j . and Garcia-Mata 
et al. [83], and in on-line notes by Caves |l6] and by Preskill |159] . 

It is prudent for students to browse among these works to find congenial points of 
view, because two of the above references are alike in the significance that they ascribe to 
the unitary invariance of quantum operations. This diversity arises because the invariance 
can be understood in multiple ways, including physically, algebraically, informatically, and 
geometrically. Our analysis will touch upon all these aspects, but much more than any of 
the above references, our approach will be geometric. 

1.3.4- Naming and applying the Theorema Dilectum It is vexing that no short name for 
the unitary invariance associated with quantum operations has been generally adopted. 
For example, this theorem is indexed by Nielsen and Chuang under the unwieldy phrase 
"theorem: unitary freedom in the operator-sum representation" |146| thm. 8.2, sec. 8.2]. 

Because we require a short descriptive name, we will call this invariance the Theorema 
Dilectum, which means "the theorem of choosing, picking out, or selecting" (from the Latin 
deligo). As our discussions will demonstrate, this name is appropriate in both its literal 
sense and in its evocation of Gauss' Theorema Egregium. 

In this article we will develop a geometric point of view in which the Theorema Dilectum 
is mainly a theorem about trajectories in state-space, and that the central practical role 
of the theorem in quantum simulations is to enable noisy quantum trajectories to be algo- 
rithmically compressed, such that efficient large-scale quantum simulation is feasible. To 
the best of our knowledge, no existing articles or textbooks have assigned to the Theorema 
Dilectum the central geometric role that this article focusses upon. 

The Theorema Dilectum is the first of two main technical terms that we will introduce 
in this review. To anticipate, the other is gabion, which is the name that we will give to 
state-space manifolds that support a certain kind of affine algebraic structure (see Figs.ls}]!] 



and Section 1.5). When gabion manifolds are endowed with a Kahler metric, we will call 
the result a gabion-Kdhler manifold (GK manifolcQ. 

GK manifolds are the state-spaces onto which we will projectively compress the quantum 
trajectories of our simulations by exploiting the Theorema Dilectum. When we further im- 
pose an antisymmetry condition upon the state-space the result is a Grassmannian gabion- 
Kdhler manifold (GGK manifold), and we will identify these manifolds as being both the 
well-known Slater determinants of quantum chemistry, and the equally well-known Grass- 
mannian varieties of algebraic geometry. 

1.3.5. Relation to geometric quantum mechanics Our recipes will embrace the strictly 
orthodox point of view that linear quantum mechanics is "the truth" to which our reduced- 
order Kahlerian state-spaces are merely a useful low-order approximation. However, at 
several points our results will be relevant to a logically conjugate point of view, known as 
geometric quantum mechanics, which is described by Ashtekar and Schilling as follows [9] 
(see also [13 172]): 

[In geometric quantum mechanics] the linear structure which is at the forefront in 
text-book treatments of quantum mechanics is, primarily, only a technical convenience 



^In a world in which every possible two-letter acronym is already in use, it is necessary to stipulate that 
this article's definition of GK manifolds does not refer to the gyrokinetic (GK) simulation codes of plasma 
physics [65] nor to the generalized Kahler (GK) manifolds of quantum field theory 
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and the essential ingredients-the manifold of states, the symplectic structure and the 
Riemannian metric-do not share this linearity. 

Thus in geometric quantum mechanics, Kahlerian geometry is regarded as a fundamental 
aspect of nature, while in our quantum MOR discussion, this same geometry is a matter 
of deliberate design, whose objective is optimizing simulation capability. Because our main 
focus is upon quantum MOR, we will comment only in passing upon those results that are 



relevant to geometric quantum mechanics {e.g., see the discussion in Section 2.12). 
1.4. Overview of the numerical simulation algorithm 

The numerical simulation algorithm of Fig. |2] is simply the formal algorithm of Fig. [T] 
expressed in a form suitable for efficient computation. Note that Fig.[2]adopts the MATLAB- 
style engineering nomenclature of model order reduction, as contrasted with the physics- 
style bra-ket notation of Fig. [T} 

The algorithm of Fig. |2] is a fairly typical example of what engineers call model order 
reduction (MOR) [3 [811 1143| \1A7\ 1148] . Rewieiiski's thesis is particularly recommended as 
a review of modern nonlinear MOR ( |164] . see also |165] ). 

1.4.1. The main ideas of projective model order reduction We will now briefly summarize 
the main ideas of projective MOR in a form that well-adapted to quantum simulation 
purposes. We consider a generic MOR problem defined by the linear equation Sip = Gtp. 
Here i/) is a state vector, 5ip is a state vector increment, and G is a (square) matrix. For 
the present it is not relevant whether 1/7 is real or complex. It commonly happens that ip 
includes many degrees of freedom that are irrelevant to the practical interests that motivate 
a simulation. 

The central physical idea of MOR is to adopt a reduced order representation ■0(c), 
where c is a vector of model coordinates, having dime ^ dimip. The central mathemat- 
ical problem of MOR is to describe the large-dimension increment 5ip by a reduced-order 
increment 6c. It is convenient to organize the partial derivatives of ip{c) as a non-square 
matrix A(c) whose elements are [A(c)]jj = dxjj^/dcj. The reduced-order increment having 
least mean-square error is obtained by the following sequence of matrix manipulations: 

(5V = GV ^ A(5c = Gt/j ^ 5c = PJ'G'ijj 5c = {A^ A)^ {K^ G ip) . (1) 

Here is the Moore-Penrose pseudo-inverse that is ubiquitous in data-fitting and model 
order reduction problems |160] , is Hermitian conjugation, and the final step relies upon 
the pseudo-inverse identity = (X^X)^X^, which is exact for any matrix X jl32| . This is 
the key step by which the master simulation equation is obtained that appears as Step B.3 
at the bottom of Fig. |2j 

The great virtue of ([l| for purposes of large-scale simulation is that (A"fA) is a low- 
dimension matrix and (AJG ip) is a low-dimension vector. Provided that both (A^^A) and 
{A!^Gip) can be evaluated efficiently, and provided also that ■j/'(c) represents the "true" 
ip with acceptable fidelity, substantial economies in simulation resources can be achieved. 
We will see that the required objectives of efficiency and fidelity both can be attained. 

1.4-2. The natural emergence of Kahlerian geometry The simulation equations, when ex- 
pressed in covariant form (Step B.3 at the bottom of Fig. |2]), provide a natural venue for 
asking fundamental geometric questions. 

For example, the low-dimension matrix d^dn = ^A^A is obviously Hermitian (whether 
-0 is real or complex) . Of what manifold is it the Hermitian metric tensor? How does this 
manifold's geometry infiuence the simulation's efficiency, fidelity, and robustness? 
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Numerical Simulation of a Quantum Trajectory 

Notation: let c„ be a vector of Kahler coordinates such 
that \ijjn) = Let be an orthonormal basis such 

that the following vectors and matrices are defined: 

(^n), = 01V'(Cn)) (H),, = {j\H\l) 

(A^,, = aO#(c))/aqlc=c„ (Mf^)).^ = OW(±)IO 
Let "define" mean "specify, but defer numeric evaluation". 
Let "increment" mean "add, but defer numeric evaluation". 
Let "compute a b" mean "evaluate a and store result in 6". 

while the simulation is running 

/* Step B.l: Specify the dynamical generator */ 
define = —iHSt/h 

/* Step B.2: Increment the dynamical generator */ 

for each of k = 1,2, ... ,m pairs of operations 

compute V'n -l^tt) ^f±)V'n ~^ Pf±)/* efficient algorithm */ 



M|V)/(pf_|_))^/^ - I [with probability p'^^^) 
Mf_)/(pf_))^/2 _ I (^otherwise) 

end for each of m pairs of operations 



increment Gn 



/* Step B.3: Increment the Kahler trajectory {ci,C2,...} */ 

compute (|A]^A^)-,7 {d(>^dK,)ji /* by an efficient algorithm */ 

compute (|A|j i/^^)^ {d(p)j /* by an efficient algorithm */ 

compute [d^dKy{d(f)) Sc^ /* is the pseudoinverse */ 

compute Cn + 5Cn Cn+l 

rescale c^+i such that xpl^_^_-^'il}n+i = 1 

end while the simulation is running 



Figure 2: Numerical algorithm for quantum model order reduction simulations. 

Steps B.1-3 are a numerical recipe that implements the simulation algorithm of Fig. [T} 
The expressions (90c?k) and {d4>) that are introduced in Step B.3 serve solely as variable 
names for the stored partial derivatives of the Kahler potential k(c,c) = ^{'il){c)\ip{c)) and 
the dynamic potential 4>{c,c) = ^{il){c)\5G\'ip{c)); it is evident that these partial derivatives 
wholly determine the simulation's geometry and dynamics. 
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To answer these questions, we will show that k is the Kdhler potential of differential 
geometry, that the metric tensor d®dK determines the Riemannian curvature of our reduced 
order state-space, and that the choice of an appropriate curvature for this state-space is 
vital to the simulation's efficiency, fidelity, and robustness. 

In preparing this article, our search of the literature did not find a previous analysis of 
MOR state-space geometry from this Riemannian/Kahlerian point of view. We did, how- 
ever, find recent work in communication theory by Cavalcante [12] and coworkers \A?>\ |2] 
that adopts a similarly geometric point of view in the design of digital signal codes. Like 
us, these authors are unaware of previous similarly geometric work ||43j "To the best of 
our knowledge this [geometric] approach was not considered previously in the context of 
designing signal sets for digital communication systems." Like us, they recognize that 
"[These state-spaces] have rich algebraic structures and geometric properties so far not 
fully explored." Also similar to us, they find fll] "The performance of a digital communi- 
cation system depends on the sectional curvature of the manifold . . . the best performance 
is achieved when the sectional curvature is constant and negative." 

Our analysis will reach similar conclusions regarding the desirable properties of nonpos- 
itive sectional curvature in the context of quantum MOR. Because the mathematical basis 
of this apparent convergence of geometric ideas between MOR theory and coding theory is 
not presently understood (by us at least), we will not comment further upon it. 

It is entirely possible that related work exists of which we are not aware. Model order 
reduction is ubiquitously practiced by essentially every discipline of mathematics, science, 
engineering, and business: the resulting literature is so vast, and the nomenclature so varied, 
that a comprehensive review is infeasible. 

It is fair to say, however, that the central role of Riemannian and Kahlerian geometry 
in model order reduction is not widely appreciated. A major goal of our article, there- 
fore, is to analyze the Riemannian/Kahlerian aspects of MOR, and especially, to link the 
Kahlerian geometry of quantum MOR to the fundamental quantum informatic invariance 
of the Theorema Dilectum. 

1.4-3. Preparing for a Kahlerian geometric analysis To prepare the way for our geometric 
analysis, at the bottom of Fig.|2]the pseudo-code defines storage variables named "(O^Sk)" 
and "((90)." For coding purposes these names are of course purely conventional (an arbitrary 
string of characters would suffice) , but these particular names are deliberately suggestive of 
partial derivatives of two scalar functions: k and (j). 

To anticipate, k will turn out to be the Kahler potential of complex differential geometry, 
which determines the differential geometry of the complex state-space, and (p will turn out 
to be a stochastic dynamical potential, which determines the drift and diffusion of quantum 
trajectories on the Kahlerian state-space. 

The link between geometry and simulation efficiency thus arises naturally because both 
the geometry and the physics of our quantum trajectory simulations are determined by the 
same two scalar functions. 

1.5. Overview of the unifying geometric ideas 

The main algebraic and geometric features of our state-space are summarized in Figs. [3]-|4j 

1.5.1. The algebraic structure of the reduced-order state space The state-space of all our 
simulations will have the algebraic structure shown in Fig. |3] We will regard this algebraic 
structure as a geometric object that is embedded in a larger Hilbert space, and we will seek 
to understand its geometric properties, including especially its curvature, in relation to our 
central topic of quantum model order reduction. 
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Figure 3: Algebraic definition of a gabion-Kaliler (GK) state-space. 

The algebraic definition of a gabion-Kahler (GK) state-space (top) ex- 
pressed equivalently as a matrix product state (MPS, bottom). By def- 
inition, the order of \'^) is the number of elements (spins) in each row's 
outer product, the rank of is the number of rows. The matrices A^™' 
that appear at bottom are, by definition, r x r matrices — hence rank 



-having diagonal elements (AW*^ 



Ikk 



and vanishing off-diagonal 



elements. Note that the matrix products are Abelian, such that the ge- 
ometric properties of the state-space are invariant under permutation of 
the spins. Note also that when the above algebraic structure is antisym- 
metrized with respect to interchange of spins (equivalent to interchange of 
columns), the state becomes a sum of Slater determinants, or equivalently 
a join of Grasssmanian manifolds (a GGK manifold). 



In the language of algebraic geometry |53t IIUU] , the geometric objects we will study are 
the algebraic manifolds that are associated with the projective algebraic varieties defined 
by the product-sums of Fig. [sj Although the literature on algebraic varieties is vast (and 
it includes many engineering applications [53j) and the literature on Riemannian sectional 
curvature is similarly vast, the intersection of these two subjects has apparently been little 
studied from an engineering point of view. This intersection, and especially its practical 
implications for quantum model order reduction, will be the main focus of our geometric 
investigations. 

The general algebraic structure of Fig.|3]is known by various names in various disciplines. 
As noted in the caption to Fig.jSj these structures are known to physicists as a matrix product 
states (often abbreviated MPS) which are widely used in condensed matter physics and ab 
initio quantum chemistry ^55] I119| I156| I157| I173| I174| I191| ; these references provide entry 
to a rapidly growing body of MPS-related literature. 
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Quantum chemists have known the algebraic structures of Fig. [3] as Hartree product 
states ^102j since 1928. Upon antisymmetrizing the outer products, we obtain the Slater 
determinants |184j that are the fundamental building-blocks of modern quantum chemistry; 
upon summing Slater determinants, and (optionally) imposing linear constraints upon these 
sums, we obtain post-Hartree-Fock quantum states \5^ . All of the theorems we derive will 



apply to Slater determinants and post-Hartree-Fock states as special cases (see Section 2.9). 
We will comment later in this section, too, upon the intimate relation of these ideas to den- 
sity functional theory (DFT). Nuclear physicists embrace these same ideas under the name 
of wave function factorization [150]. Beylkin and Mohlenkamp [TO note that statisticians 
call essentially the same mathematical objects canonical decompositions and also parallel 
factors. 

As Leggett and co-authors have remarked with regard to the similarly immense literature 
on two-state quantum systems: "The topic of [this] paper is of course formally a problem 
in applied mathematics. . . . Ideas well known in one context have been discovered afresh in 
another, often in a language sufficiently different that it is not altogether trivial to make the 
connection. ... [In such circumstances] the primary purposes of citations are to help the 
reader understand the paper, and the references in the text are chosen with this in mind" 
|129] . These same considerations will guide our discussion. 

The general utility of affine algebraic structures for modeling purposes first came to 
our attention in a highly readable Mathematical Intelligencer article by Mohlenkamp and 
Monzon j l40| ; two further articles by Beylkin and Mohlenkamp [15', TB^ are particularly rec- 
ommended also. Beylkin and Mohlenkamp call the algebraic structure of Fig. |3] a separated 
representation, and they have this to say about it |15] : 

When an algorithm in dimension one is extended to dimension d, in nearly every case 
its computational cost is taken to the power d. This fundamental difficulty is the 
single greatest impediment to solving many important problems and has been dubbed 
the curse of dimensionality. For numerical analysis in dimension d, we propose to 
use a representation for vectors and matrices that generalizes separation of variables 
while allowing controlled accuracy. . . . The contribution of this paper is twofold. First, 
we present a computational paradigm. With hindsight it is very natural, but this 
perspective was the most difficult part to achieve, and it has far-reaching consequences. 
Second, we start the development of a theory that demonstrates that separation ranks 
are low for many problems of interest. 

In a subsequent article Beylkin and Mohlenkamp go on to say [T6^ "The representation 
seems rather simple and familiar, but it actually has a surprisingly rich structure and is 
not well understood." These remarks are remarkably similar in spirit to the coding theory 
observations of Cavalcante et al. that were reviewed in Section II. 4[ 

For us, Kahlerian algebraic geometry will provide a shared foundation for understanding 
the accelerating progress that all of the above large-scale computational disciplines have 
witnessed in recent decades. 

1.5.2. The medieval idea of a gabion, and its mathematical parallels Deciding what to call 
the geometric state-space of quantum model order reduction is a vexing problem. We have 
seen that various plausible names include "Hartree products," "Slater determinants," "sep- 
arated representations," "matrix product states," "wave function factorizations," "product- 
sum states," "canonical decompositions," and "parallel factors." 

A shared disadvantage of the above names is that there is no precedent for associating 
them with the geometric properties that are the main focus of our investigations. We 
therefore seek an encompassing name for these state-spaces viewed as geometric entities. 
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Finding no precedent in the literature, and desiring a short name having a long-estabhshed 
etymology, we will call them by the medieval name of gabions |192] , or more formally, gabion 
manifolds (this name arose spontaneously in the course of a seminar). 

Most readers will have seen gabions numerous times, perhaps without recognizing that 
they have a well-established name. "Gabion" is the generic engineering name for a mesh 
basket that is filled with a weighty but irregularly-shaped material such as rocks or lumber, 
then stacked for purposes of reinforcement, erosion control, and fortification. In medieval 
times gabions were made of wicker or reed; Fig. QA) shows a typical medieval gabion. We 
will see that the defining geometry property of gabion manifolds is that they are possessed 
of a web of geodesic lines that constrain the curvature of the manifold, rather as the wicker 
reeds of a physical gabion constrain the curved rocks and boulders held inside. Like physical 
gabions, gabion manifolds come in a wide variety of sizes and shapes that are suitable for 
numerous practical purposes. 

We will postpone giving a formal — and necessarily rather abstract — mathematical def- 



inition of a gabion until Section 2.5 For the present our main objective is to informally 
describe the geometric properties that will motivate this formal definition. 

1.5.3. The geometric properties of gabion-Kdhler (GK) manifolds The main geometric 
properties of GK manifolds that are relevant to quantum simulation are depicted in Fig.|4]^B- 
H). We will now survey these properties, and in doing so, we will introduce some of the 
nomenclature of Kahlerian geometry. 

We begin our geometric overview by remarking that even though Hilbert space is a com- 
plex state-space, a common viewpoint among mathematicians is that a complex manifold 
is a real manifold that is endowed with an extra symmetry, called its complex structure (see 



Sections 2.2 and 2.6 for details). For purposes of our geometric analysis, we will simply 



ignore this complex structure until we are ready to apply the quantum Theorema Dilectum 



(Section 2.6). Until then we will treat gabion manifolds as real manifolds. 

In particular, the state-space of quantum mechanics has a natural real-valued measure 
of length. Specifically, along a time-dependent quantum trajectory |^(t)) it is natural to 
define a real- valued velocity v{t) whose formal expression can be written equivalently in 
several notations: 

dimH/2 

7;(i)2 = g(V'(t),#)) = (V'(t)|#)> = -0(O>W= Yl (2) 

1=1 

Here ip{t) = dip{t)/dt, and we have used first the abstract notation of differential geome- 
try (in which g{. ■ ■) is a metric function), then the Dirac bra-ket notation of physics, then 
the matrix-vector notation of engineering and numerical computation, and finally the cum- 
bersome but universal notation of components and sums over indices. We will assume an 
entry-level familiarity with all four notations, since this is a prerequisite for reading the 
literature. 

As a token of considerations to come, the factor of dim 7^/2 in the index limit of ^ 
above arises because we will regard a complex manifold like as being a real manifold 
of dimension 2n. Thus we will regard the complex plane <L as a two-dimensional (real) 
manifold, and an spin-1/2 quantum state as a point in a Hilbert space TC having diniTC = 4 
(real) dimensions. This viewpoint leads to an ensemble of conventions that we will review 



in detail in Section 2.6 For now, we note that the arc length s along a trajectory is 
s = J v{t) dt, so that geometric lengths in quantum state-spaces are dimensionless. An 
equivalent differential definition is to assign a length increment ds to a state increment \dtp) 
via {dsf = (diljldip). 
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(g) robust projective fidelity (h) tlie Theorema Egregium (i) the Theorema Dilectum 




Figure 4: Geometric principles of quantum model order reduction (QMOR). 
See Section 1.5 for a discussion of these principles. 
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Since we can now compute the real-valued length of an arbitrary curve on the gabion 
manifold, all of the usual techniques of differential geometry can be applied, without special 
regard for the fact that the state-space is complex. 

1.5.4- GK manifolds are endowed with rule fields Beginning a pictoral summary of our 
geometric results, we first note that state-space gabions resemble physical gabions in that 
they are naturally endowed with a geometric mesh, which is comprised of a network of lines 
called rules, as depicted in Fig. |4]^B). More formally, they are equipped with vector fields 



having certain mathematical properties (see (14) and Definition 2.1) such that the integral 



curves of the rule fields have the depicted properties. 

Postponing a more rigorous and general definition of gabion manifolds until later (see 



Section 2.5 ), we can informally define a gabion rule to be the quantum trajectory associated 
with the variation of a single coordinate "c^ in the algebraic structure of Fig. [sj holding 
all the other coordinates fixed to some arbitrary set of initial values. We see that gabion 
rules are rays (straight lines) in the embedding Hilbert space, and hence, the gabion rules 
are geodesies (shortest paths) on the gabion manifold itself. As depicted in Fig. |4]^C), the 
set of all gabion points that belong to a rule is (trivially) the set of all gabion points itself, 
which is the defining characteristic of a gabion being ruled. Furthermore, we will show that 
at any given point, the vectors tangent to the rules that pass through that point are a basis 
set. 

1.5.5. GK geometry has singularities Are gabion manifolds geometrically smooth, or do 
they have singularities? As depicted in Fig. |4|^D), we will show that gabion manifolds have 
pinch-like geometric singularities. Algebraically these singularities appear whenever two 
or more rows of the product-sum in Fig. [3] are equal. Geometrically, we will show that 
the Riemann curvature diverges in the neighborhood of these singular points. However, 
it will turn out that the continuity of the geodesic rules is respected even at the singular 
points, so that gabion manifolds are geodesically complete. Pragmatically, this means that 
our numerical simulations will not become "stuck" at geometric singularities. 

1.5.6. GK projection yields compressed representations As depicted in Fig. |4]^E) model 
order reduction is achieved by the high-fidelity projection of an "exact" state in the 
large-dimension Hilbert space onto a nearby point \ip]c) of the small-dimension gabion. It 
can be helpful to view this projection as a data-compression process. By analogy, the state 
1^) is like an image in TIFF format; this format can store an arbitrary image with perfect 
fidelity, but consumes an inconveniently large amount of storage space. The projected state 
IV'a:) on the gabion is like an image in JPEG format; lesser fidelity, but good enough for 
many practical purposes, and small enough for convenient storage and manipulation. We 
thus appreciate that data compression can be regarded as a kind of model order reduction; 
the two processes are fundamentally the same. 

1.5.7. GK manifolds have negative sectional curvature Practical quantum simulations re- 
quire that the computation of order- reducing projections be efficient and robust, just 
as we require image compression programs to be efficient and robust. As depicted in 
Fig. |4]^F), order-reduction projection becomes ill-conditioned when the state-space man- 
ifold is "bumpy" , in which case a numerical search for a high-fidelity projection can become 
stuck at local minima that yield poor fidelity. We will prove that the presence of a ruled 
net guarantees that gabion manifolds are always smooth rather than bumpy. 

Resorting to slightly technical language to say exactly what we mean when we assert 
that gabion manifolds are not bumpy, in our Theorem |2.1| we will prove that a gabion 
has nonpositive sectional curvature for all sections on its geodesic net. This means that 
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gabion manifolds can be envisioned as a net of surfaces that have the special property of 
being saddle-shaped everywhere (as contrasted with generic surfaces having dome-shaped 
"bumps"). As depicted in Fig. |4]^G), the saddle-shaped curvature helps ensure that order- 
reducing projection onto gabion manifolds is a numerically well-conditioned operation. 

1.5.8. GK manifolds have an efflorescing global geometry As depicted in Fig. |4|h), when 
the number of state-space dimensions becomes very large, it becomes helpful to envision 
nonpositively curved manifolds as flower-shaped objects composed of a large number of 
locally Euclidean "petals." This physical picture has been vividly conveyed by recent col- 
laborative work between mathematicians and fabric artists [12]; the work of Taimina and 
Henderson on hyperbolic manifolds is particularly recommended [ 109j . 

When working in large-dimension spaces we will heed also Dantzig's remark that "one's 
intuition in higher dimensional space is not worth a damn!" [5 . For purposes of quantitative 
analysis we will rely upon Gauss' Theorema Egregium |86J to analyze the Riemannian 
and Kahlerian geometric properties of gabion manifolds. We will prove that "number of 
petals" becomes exponentially large, relative to dim/C, such that the petals loosely fill the 
embedding Hilbert space. In this respect our geometric analysis will parallel the informatic 
analysis of Nielsen and Chuang |146j : their Fig. 4.18 is broadly equivalent to our Fig. |4|^H). 

Our analysis will therefore establish two geometric properties of gabion manifolds: they 
are strongly curved, and they are richly endowed with straight-line rules. We will show that 
these gabion properties are essential to the efficiency, robustness, and fidelity of large-scale 
MOR. Later on, in Sections |4. 6. 4f|4".6.6[ we will establish a relation between these properties 
and compressive sampling (CS) theory. 

1.5.9. GK basis vectors are over-complete We will nowhere assume that the basis vectors of 
the underlying algebraic structure of Fig. [3] are orthonormal; they might refer for example 
to the non-orthonormal gaussian basis states of quantum chemistry. A major geometric 
theme of our analysis, therefore, is that the negative sectional curvature of gabion manifolds 
helps generically account for the observed efficiency, fidelity, and robustness of gabion-based 
modeling techniques in many branches of science, engineering, and mathematics. 

1.5.10. GK manifolds allow efficient algebraic computations Upon restricting our attention 
to the special case of gabion-Kahler (GK) manifolds, we will show that the existence of a 
ruled geodesic net allows the sectional curvature and the Riemann curvature tensors of GK 
manifolds to be calculated easily and efficiently. To anticipate, we will present data from 
Riemann curvature tensors having dimension up to 188, which we believe are the largest- 
dimension curvature tensors yet numerically computed. We will see that it is Kraus' "long 
list of miracles" that makes large-scale numerical curvature computations feasible, and that 
these same miracles are equally essential to large-scale quantum dynamical calculations. 

1.5.11. GK manifolds support the Theorema Dilectum One geometric idea remains that is 
key to our simulation recipes. For gabion manifolds to represent quantum trajectories with 
good fidelity, some physical mechanism must be invoked to compress quantum trajectories 
onto the petals of the gabion state-space. That key mechanism is, of course, the Theorema 
Dilectum that was mentioned in Section [l.3.4[ 

From a geometric perspective, the Theorema Dilectum guarantees that noise can always 
be modeled as a measurement process that acts to compress trajectories onto the GK 
petals. As depicted in Fig. |4]^I), quantum simulation can be envisioned geometrically as a 
process in which compression toward the GK petals, induced by measurement processes, 
competes with expansion away from the petals, induced by quantum dynamical processes. 
The balance of these two competing mechanisms determines the MOR dimensionality that is 
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required for good fidelity — the "petalthickness." Algebraically this petal thickness increases 
in proportion to the rank of the product-sum algebraic structure of Fig. [3} 

Thus for us, trajectory compression is not a mathematical "trick," but rather is a rea- 
sonably well-understood and well-validated quantum physical mechanism, originating in the 
Theorema Dilectum, that compresses quantum trajectories to within an exponentially small 
fraction of the Hilbert phase space. 

This noise-induced trajectory compression is the loophole by which QMOR simulations 
evade the no-simulation arguments of Feynman [67^, as reviewed by Nielsen and Chuang 
|146| see their Section 4.7]. 

1.5.12. GK manifolds support thermal equilibria We will see that this covert-measurement 
approach encompasses numerical searches for ground states. Specifically, by explicit con- 
struction, we will show that contact with a zero-temperature thermal reservoir can be 
modeled as an equivalent process of covert measurement and control, in which the role of 
"temperature" is played by the control gain, such that zero temperature is associated with 
optimal control. 

From this QMOR point of view, the calculation of a ground-state quantum wave func- 
tion is a special kind of noisy quantum simulation, in which noise is present but masked 
by optimal control. This is how QMOR reconciles the strong arguments for the general 
infeasibility of ah initio condensed-matter calculations (as reviewed by, e.g.^ Kohn |121] ) 
with the widespread experience that numerically computing the ground states of condensed 
matter systems is often, in practice, reasonably tractable |80j . 

1.5.13. GK manifolds support fermionic states Readers familiar with ab initio quantum 
chemistry, and in particular with density functional theory (DFT) |28| I114| I12H see Cappele 
[39j for an introduction] will by now recognize that QMOR and DFT are conceptually parallel 
in numerous fundamental respects: the central role of the low-dimension Kahler manifold 
of QMOR parallels the central role of the low-dimension density functional of DFT; the 
closed-loop measurement and control processes of QMOR parallel the iterative calculation 
of the DFT ground state; QMOR's fundamental limitation of being formally applicable 
only to noisy quantum systems parallels DFT's fundamental limitation of being formally 
applicable only to ground states; QMOR and DFT share a favorable computational scaling 
with system size. 

Yet to the best of our knowledge — and surprisingly — the geometric techniques that that 
this article will deploy in service of QMOR have not yet been applied to DFT and related 
techniques of quantum chemistry and condensed matter physics [80] . A plausible starting 
point is to impose an antisymmetrizing Slater determinant-type structure upon the algebraic 
outer products of Q. 

Some analytic results that we have obtained regarding the Kahlerian geometry of Slater 
determinants are summarized in Section 2.9, With further work along these lines, we believe 
that there are reasonable prospects of establishing a geometric/info rmatic interpretation, 
via the Theorema Egregium and the Theorema Dilectum, of the celebrated Hohenberg- 
Kohn and Kohn-Sham Theorems of DFT |114j and their time-dependent generalization the 
Runge- Gross Theorem |28j . 

A physical motivation for this line of research is that the Theorema Dilectum of QMOR 
and the Hohenberg-Kohn Theorem of DFT embody essentially the same physical insight: the 
details of exponentially complicated details of quantum wave functions are only marginally 
relevant to the practical simulation of both noisy systems (QMOR) and systems near their 
ground-state (DFT). 

At present, the two formalisms differ mainly in their domain of application: QMOR is 
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well-suited to simulating spatially localized systems at high temperature {e.g., spin systems) 
while DFT is particularly well-suited to simulating spatially delocalized systems {e.g., mol- 
ecules and conduction bands) at low temperature. In the future, as QMOR is extended 
to delocalized systems while the methods of quantum chemistry are increasingly extended 
to dynamical systems [28! 180], opportunities will in our view arise for cross-fertilization 
of these two fields, both in terms of fundamental mathematics and in terms of practical 
applications. 

1.6. Overview of contrasts between quantum and classical simulation 

In aggregate, the formal, numerical, algebraic, and geometric concepts summarized in the 
preceding sections and in Figs. [lj|4] are in many respects strikingly parallel to similar con- 
cepts in the computational fluid dynamics (CFD), solid mechanics, combustion theory, and 
many other engineering disciplines that entail large-scale simulation using MOR. 

However, it is evident that quantum MOR is distinguished from real-valued (classical) 
MOR by at least four major differences, which we will now summarize. 

1.6.1. The Theorema Dilectum is fundamental and universal The first difference is that 
the Theorema Dilectum describes an invariance of quantum dynamics that is fundamental 
and universal. Its physical meaning, as we will see, is that it enforces causality. Nonlinear 
classical system do not possess any similarly universal invariance, which is in our view 
a major contributing reason that "developing effective and efficient MOR strategies for 
nonlinear systems remains a challenging and relatively open problem" |164| p. 20]. 

Our results, both analytical and numerical, will suggest that noisy quantum systems 
are fundamentally no harder to simulate than nonlinear classical systems, provided that 
the Theorema Dilectum is exploited to allow high-fidelity dynamical projection of quantum 
trajectories onto a reduced-order state-space. 

1.6.2. Quantum state-spaces are veiled The second difference is a consequence of the first. 



As discussed in Sec. 1.3 to fully exploit the power of the Theorema Dilectum we are re- 
quired to embrace the ultra-orthodox principle of never looking at the quantum state space. 
Furthermore, when we examine classical state-spaces more closely, we find that they too 
are encumbered with ontological ambiguities that precisely mirror the "spooky mysteries" 



of quantum state-spaces. As discussed in Sections |3. 2. 7[ this modern recognition of spooky 
mysteries in classical physics echoes work in the 1940s by Wheeler and Feynman |1961I197] . 

1.6.3. Noise makes quantum simulation easier The third difference is that higher noise lev- 
els are beneficial to QMOR simulations, because they ensure stronger compression onto the 
GK petals, which allows lower-rank, faster-running GK state-spaces to be adopted. Later 
we will discuss the interesting question of whether this principle, together with the con- 
comitant principle "never look directly at the quantum state-space," have classical analogs. 
We will tentatively conclude that the Theorema Dilectum does have classical analogs, but 
that the power of this theorem is much greater in quantum simulations than in classical 
ones. 

1.6.4. Kdhlerian manifolds are geometrically special Broadly speaking, Kahlerian geome- 
try is to Riemannian geometry what analytic functions are to ordinary functions. This addi- 
tional structure is one of the reasons why the mathematician Shing-Tung Yau has expressed 
the view |200| p. 46] "The most interesting geometric structure is the Kahler structure." 
From this point of view, the geometry of real-valued MOR state-spaces is mathematically 
interesting, and the analytic extension of this geometry to Kahlerian MOR state-spaces is 
even more interesting. 
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Let us state explicitly some of the analogies between analytic functions and Kahler man- 
ifolds. We recall that generically speaking, analytic functions have cuts and poles. These 
cuts and poles are of course exceedingly useful to scientists and engineers, since they can 
be intimately linked to physical properties of modeled systems. Similarly, the GK mani- 
folds that concern us have singularities, as depicted in Fig. |4]^D). Physically speaking, they 
are associated with regions of quantum state-space that locally are more nearly "classical" 
than the surrounding regions, in the sense that the local tangent vectors that generate high- 
order quantum correlations become degenerate. It is fair to say, however, that the deeper 
geometric significance of Kahlerian MOR singularities remains to be elucidated. 

Just as contour integrals of analytic functions can be geometrically adjusted to make 
practical reckoning easier, we will see the Theorema Egregium allows the trajectories arising 
from the drift and diffusion of noise and measurement models to be geometrically (and 
informatically) adjusted to match state-space geometry, and thereby improve simulation 
fidelity, efficiency, and robustness. 

More broadly, Yau notes [200. p. 21]: "While we see great accomplishments for Kahler 
manifolds with positive curvature, very little is known for Kahler manifolds [having] strongly 
negative curvature." It is precisely these negatively-curved Kahler manifolds that will con- 
cern us in this article, and we believe that their negative curvature is intimately linked 
to the presence of the singularities mentioned in the preceding paragraph. We hope that 
further mathematical research will help us understand these connections better. 

2. The sectional curvature of gabion-Kahler (GK) state-spaces 

We will now proceed with a detailed derivation and analysis of our quantum simulation 
recipes. Our analysis will "unwind" the preceding overview: first we analyze the geometry 
of Fig. |4j as embodying the algebraic structure of Fig. [3j using the numerical techniques of 
Fig. [2] Only at the very end will we calibrate our recipes in physical terms, via the quantum 
physics of Fig. [T} 

2.1. Quantum MOR state-spaces viewed as manifolds 

To construct our initial example of a gabion state-space, we will consider the following alge- 
braic function ■j/'(c), whose domain is a four-dimensional manifold of complex coordinates 
c = {c^, c^, , c^} and whose range in a four-dimensional Hilbert space is the set of points 
that can be algebraically represented as {c^c", c^c^, c^c^, c^c^}. In the notation of Fig. [2] 
this function is 



^(c) 



V" 




C3" 






c\ 



' ii\ - c^c^ = 

- c^d^ = 

V'3 - C^C' = 

[ V4 - c^c^ = 



2„3 _ n ' (^) 



where "0" is the outer product. The superscripts on the c* variables are indices rather 
than powers, as will be true throughout this section. From an algebraic geometry point 
of view, ^ defines a projective algebraic variety [53j (also called a homogeneous algebraic 
variety) over variables {ipi. z G 1,4} that is specified above in parametric form in terms of 
parameters {c*: i G 1,4}. 

By definition, our example of a gabion state-space manifold is the solution set of this 
algebraic variety, and thus our state space is an algebraic manifold. Physically speaking, 
i/^(c) is the most general (unnormalized) quantum state of two spin 1/2 particles sharing 
no quantum entanglement. 

We will now show that this state-space is a Kahlerian manifold that has negative sec- 
tional curvature (under circumstances that we will describe) and that this property is ben- 



2.2 Regarding gabion manifolds as real manifolds 
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eficial for simulation purposes (for reasons that we will describe). 

Practical computational considerations: We begin by remarking that the basic algebraic 
construct "argl (E) arg2" that appears in ([3| can be readily implemented by the built-in 
functions of most scientific programming languages and libraries; for example in MATLAB 
by the construct "reshape ( (argl*arg2 ' ) ' , [] ,1)" and in Mathematica by the construct 
"Outer [Times , argl , arg2] //Flatten" . 

Similar idioms exist for the efficient evaluation of more complex product-sum structures. 
Although we will not describe our computational codes in detail, they are implemented in 
MATLAB and Mathematica in accord with the general ideas and principles for efficient 
addition, inner products, and matrix-vector multiplication that are described by Beylkin 
and Mohlenkamp [16 . 

The abstract geometric point of view: From an abstract point of view, the algebraic structure 
Q can be regarded as a sequence of maps 

_ surjective ,^ injective ^ , , , > 

c /c ^ n, (4) 

where C is the manifold of complex variables {0^,0^,0^,0^}, the gabion manifold fC is the 
range of ■0 in TC, and TC is the larger Hilbert space within which IC is embedded. To 
appreciate the surjective and injective nature of these maps, we notice in Q that i/'(c) is 
invariant under {0^,0^,0^,0^} — > {1, c^/c^, c^c^, c^c^}. More generally, it is clear that one 
coordinate can be set to any fixed nonzero value, without altering ip{c), by an appropriate 
rescaling of the other three variables. 

In our example ([3|, the dimensions of the three manifolds C, JC, and 7i are therefore 

dimC = 2x4 = 8, dim/C = 2x3 = 6, dimW = 2 x 4 = 8, (5) 

where the factors of two arise because these are complex manifolds. We see that the map 
C — > /C is surjective (because dim/C < dimC), while /C — > W is injective (because IC is 
immersed in TC). 

2.1.1. Defining gabion pseudo-coordinates We will call the variables {0^,0^, c^, c^} pseudo- 
coordinates. They are not ordinary coordinates because C — > /C is surjective rather than 
bijective, or to say it another way, open sets on C are not charts on IC. Whenever we require 
an explicit coordinate basis, we can simply designate any one c'^ to be some arbitrary fixed 
(nonzero) value, and take the remaining {c* : i 7^ k} to be coordinate functions. 

In practical numerical calculations — where these algebraic structures are called "sepa- 
rated representations," "matrix product states," or "Slater determinants" — pseudocoordinate 
representations are adopted almost universally. Therefore, we will sometimes simply call 
the c's "coordinates" ; this will make it easier to link the numerical algorithm of Fig. [2] to 
the geometric properties of IC. 

2.2. Regarding gabion manifolds as real manifolds 

Now we will begin analyzing in detail the curvature of the gabion manifold fC. For geometric 
purposes it is convenient to regard Ti not as a complex vector space, but as a Euclidean 
space, such that i/> is a vector of real numbers that in our simple example has the eight 
components {V''"} = {3?(^/'^), . ■ ■ , ■ . . , 9(V'^)}. Similarly, we specify real co- 

ordinates on C via 0^ = -\- iy^, and with a long-term view toward interfacing with the 
Kahler geometry literature we agree to specify these real coordinates in the conventional 
order {x^ , . . . ,x^,y^ , . . . ,y'^} = {r^,...,r^} = {r"-}. Thus {d/dr"-} is a complete set of 
vectors on C. 



2.3 "Push-hutton" strategies for curvature analysis 
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2.2.1. Constructing the metric tensor Then the map tp : C ^ 7i induces a metric tensor g 
upon C via the Euchdean metric of 7i. The components of g evidently are 



9ab = 9 



fd d \ 


5'0(c(r)) 




3'j/'(c(r)) 




dr"- 







(6) 



This also suffices to define g as a metric tensor on /C, provided we restrict our attention — as 
for MOR purposes we always will — to functions on /C having the functional form /(i/'(c(r))), 
such that the tangent vectors {d/dr'^} always act either directly or indirectly (via the chain 
rule) upon '0(c(r)). Then knowledge of g allows us to compute via ^ the velocities and 
path lengths of arbitrary trajectories on /C, as is required of a metric on /C. 

2.2.2. Raising and lowering the indices of a pseudo- coordinate basis Considered as a co- 
variant matrix, the indices of gab range over an over-complete basis set, and therefore gab 
is singular. It follows that we cannot construct a contravariant matrix g""^ in the usual 
manner, by taking a matrix inverse of gab- To evade this difficulty we define the contravari- 
ant metric tensor to have components g"'^ = {QabY ■, where is the same Moore-Penrose 
matrix pseudoinverse that appears in Step B.3 of Fig. [2] and that was discussed following 
Eq. Q. 

It is easy to verify that g""^ and Qab act to raise, lower, and contract tensor indices in 
the usual manner, with a single important difference: the operation of raising followed by 
lowering is no longer the identity operator, but rather is a projection operator, in conse- 
quence of the general pseudoinverse identity {X^X)"^ = X''X. Physically this projection 
annihilates tangent vectors on /C whose length is zero. 

2.2.3. Constructing projection operators in the tangent space From these identities it fol- 
lows that at a specified point tp^^ of /C, the local operator Pic{tpic) that projects vectors in 
Ti. onto the tangent space of /C at \ip) is 



K) 



dimC 

E 

a,6=l 



^ab 



Q^b 



(7) 



In the interest of compactness, we will often write Pjq rather than Pn^^ipj^). We readily 
verify that the projective property PjcPfc = Pjc follows from the definition of g"''' given in 
Q and the general pseudo-inverse identity X^XX^ = X^. 

The ability to construct the projection P/c solely from tangent vectors and the local 
metric tensor g will play a central role in our geometric analysis of /C. 

2.3. "Push-button" strategies for curvature analysis 

At this point we can analyze /C's intrinsic geometry by either of two strategies. The first 
strategy, which can be wholly automated, is to fix in our example problem (say) r"^ = 1 
and = so that the remaining {r^,r^, . . . ,r^} can be regarded as conventional coordi- 
nate functions on the six-dimensional gabion /C. The now-restricted set of tangent vectors 
associated with {r^, r^, . . . , r^} constitutes a coordinate basis, such that ^ specifies the 
metric tensor for this basis. By construction, this metric has no null vectors, and hence gab 
is invertible. 

The intrinsic geometric properties of IC can then be automatically computed by any of 
the many symbolic manipulation packages that are available for research in general rela- 
tivity. This automated approach allows us to "push the button" and discover that for our 
simple example the scalar Riemann curvature of /C is given by the remarkably simple 
expression R = —8/{tp-ip). 



2.4 The sectional curvature of gabion state-spaces 
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2.3.1. The deficiencies of push-button curvature analysis What is unsatisfying about an 
automated coordinate-based analysis, however, is that this simple integer result is obtained 
as the result of seemingly miraculous cancellations of high-order polynomials. This method 
produces no insight as to why such a simple integer result is obtained, or why the sign of 
the curvature is negative, or whether this simplicity is linked to /C's ruled structure. 

Another objection to coordinate-based analysis is that it forces us to "break the sym- 
metry" of the coordinate manifold C by designating arbitrary fixed values for arbitrarily 
selected gabion coordinates. This is undesirable because our quantum simulation algo- 
rithms in Figs. [T}|2] do not break this symmetry. To do so in our geometric analysis would 
unnecessarily obstruct our goal of linking quantum simulation physics to the Kahlerian 
geometry of /C. 

We will therefore develop a Riemannian/Kahlerian curvature analysis of the gabion 
manifold /C that fully respects the algebraic symmetries, not of /C, but of C. For this 
purpose the sectional curvature proves to be an ideal mathematical tool. 

2.4- The sectional curvature of gabion state-spaces 

Because /C has a natural embedding in the Euclidean manifold TC, our analysis of sectional 
curvature is able to follow quite closely the embedded geometric reasoning of Gauss' original 
derivation of the Theorema Egregium [86]. This approach has the advantage of yielding 
immediate physical insight. Equally important, this approach can be readily adapted to 
accommodate the pseudo-coordinate tangent basis that is most natural for analyzing gabion 
geometry. 

As depicted in Fig. |4]^E), we choose an arbitrary point on K, and define tangent vectors 
U and V on fC to be directional derivatives 

dim C r\ dim C r\ 

U-y^u""^ and V=y^v''^. (8) 

a=l a=l 

Because the map C ^ /C is surjective, the representation of U and y as a sum over compo- 
nents and v"" is nonunique, and we will take care to establish that our sectional curvature 
calculations are not thereby affected. 

2.4-1- Remarks on gabion normal vectors Conjugate to the tangent space of /C at a given 
point is the space of vectors in Ti. that are normal to the tangent space. We specify n to be 
an (arbitrarily chosen) unit vector in that normal space, i.e., to be a vector satisfying 

dip(r) 

h ■ h = 1 and h • ^ ^ = n ■ ip^^ = , (9) 

where we have adopted the usual notation that a comma preceding a subscript (s) indicates 
partial differentiation with respect to the indexed variable(s). The sign of h will not be 
relevant. We remark that h is not unique because the codimension of /C (by definition 
codim/C = dim'H — dim/C) is in general greater than unity. Looking ahead, in some of our 
large-scale numerical examples codim/C will be very large indeed, of order 2 x 2^^ ~ 512, 000. 

2.4-2. Computing the directed sectional curvature With reference to the vectors U, V, and 
h depicted in Fig.|4](E), we define a scalar function S{U, V, h), which we will call the directed 
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sectional curvature, to be 



S{U,V,n) 



dimC 

E 

a,b,c,d=l 



n ■ 1p,ab 
n ■ il^.bd 



dimC 

E 

a,b,c,d=l 



(10) 



yaybyCyd 



Here |. . .| denotes the determinant. 

If we recall that h ■ rp^^ = 0, per ([o]), then it is straightforward to verify that 5 is a scalar 
under coordinate transformations. It further satisfies the identity 



S{U, V, h) = S{aU + 13V, jU + 6V, h) 



(11) 



for arbitrary real a, /?, 7, and 6. Thus S{U, V, fi) is a real-valued geometric invariant of the 
two-dimensional tangent subspace spanned by U and V . 

In the preceding paragraph we emphasize that a, /?, 7, and 5 are real-valued because 
later on, when we admit a complex structure, it will not be true that phase-shifting U 



and/or V leaves the sectional curvature invariant (see Section 2.7). 



The denominator of (10) has a simple physical interpretation as the geometric area of 
the section defined by U and V] this quantity is often written as |C/ A V\'^. In terms of the 
metric function g{U, V) = i,gabu"v^ we have 



\U ^V\^ = g{U,U)g{V,V) - g{U,Vf 



(12) 



2.4-3. Physical interpretation of the directed sectional curvature For a two-dimensional 
surface embedded in a three dimensional space — the case considered by Gauss — the above 
expression reduces to the familiar expression S{U, V,h) = l/(i?ii?2), where Ri and R2 are 
the principal radii of curvature of the surface. In higher dimensions S{U, V, h) describes 
the Gaussian curvature of a two-dimensional section of /C — a two-dimensional submanifold 
that is locally tangent to U and V — that has been projected onto the three-space spanned 
hy {U,V,n}. 

For MOR purposes, this means that whenever S{U, V, h) is negative we are guaranteed 
local concavity of the state-space /C as viewed along n {i.e, as viewed "from above") along at 
least one curve that is locally tangent to some linear combination of U and V. The resulting 
physical picture is shown in Fig. |4]|^G). For MOR purposes our goal will be, therefore, 
to choose state-space manifolds such that S{U, V, h) is negative, in the expectation that 
the associated local concavity of /C will improve the robustness of projective model order 
reduction. 

2.4-4- Definition of the intrinsic sectional curvature Mathematicians usually prefer to de- 
scribe the curvature of IC in intrinsic terms. To accomplish this it is convenient to sum 
over a complete orthonormal set {ui} of vectors tangent to IC. We use the identity 
^ ■ n (g) n = Pic = I — P/c, where / is the identity operator and P/c is the projection 
matrix given in ([7|, to obtain 



S{U,V) 



CO dim JC 

S{U,V,n. 

i=l 
dimC 



E 



\U AV\^ 



cd 



^a^b^c^d 



(13) 
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Now all reference to unit normals has disappeared, because we have already established in 
that P,c 

can be described in intrinsic terms. However, the above expression still refers to 



the embedding Hilbert space via ip. It will not be until later on (specifically, following (26)) 
that we show that S(U, V) is determined solely by the metric tensor and its derivatives. 

2.5. The formal definition of a gabion manifold 

For general tangent vectors U and V, the directed sectional curvature S{U, V, h) and the 
intrinsic sectional curvature S{U, V) can be either positive or negative. We will now derive 
a condition on U and V which is sufficient for S{U, V, n) and S{U, V) to be nonpositive, 
and we will use this condition to motivate a formal definition of a gabion manifold. 

We recall that /C is a reduced-dimension state-space manifold that is embedded in a 
larger- dimension Euclidean manifold Ti.. Thus each individual component of the state- 
vector ■0 of defines a scalar function on /C. If we wish, we may regard the metric g of (§ 
and the normal vector n of ([9| as intrinsically defined in terms of the scalar functions i/'; 
this eliminates any formal reference to the embedding manifold 7i. We define a rule vector 
field, or simply rule field, to be any vector field V on /C satisfying 

VyVv^ = 0. (14) 

The motivation for this definition is simply that the above equation is both intrinsic and 
geometrically covariant, and furthermore, it is manifestly satisfied by the vector field that 
is associated with each gabion pseudo-coordinate; these coordinates thus are canonical 
examples of rule fields. We define rule lines, or simply rules, to be the integral curves of a 
rule field. It is straightforward to show that rule lines are geodesies; this formally justifies 
our earlier depiction of rules as "straight lines" in Fig. |4] and in Section 1.5, We define rule 



tangent vectors, or simply rule vectors, to be vectors that are locally tangent to a rule line. 
We now formally define a gabion manifold as follows: 

Definition 2.1. A gabion manifold is a manifold endowed with rule fields whose rule tangent 
vectors constitute a local basis at every point of the manifold. 

This definition is more restrictive than the usual definition of a ruled manifold [7T], in 
which there is no requirement that the rule vectors provide a local basis. Roughly speaking, 
therefore, a gabion manifold is a ruled manifold that is exceptionally rich in rule structure. 

Associated with a rule vector V, we define local rule coordinates such that VyVyi/? = 
takes the component form 'Yliab'^"^'^^'^ ,ab — 0- Thus gabion pseudo-coordinates are local 



rule coordinates. Evaluating (13) in local rule coordinates, we see that whenever either U 
or y is a rule vector, the first numerator term vanishes, and the remaining numerator term 
is nonpositive. This proves the following theorem: 

Theorem 2.1. Let U be a rule vector at an arbitrary point on a manifold fC, let V be an 
arbitrary tangent vector at that same point, and let h be an arbitrary unit vector normal 
to the tangent space. Then the directed sectional curvature satisfies S{U, V,n) < and 
therefore, the intrinsic sectional curvature satisfies S{U,V) < 0. 

In physical terms, any two-dimensional section of /C that includes a rule vector has neg- 
ative sectional curvature. Since for gabion manifolds, the local rule tangents form an 
over-complete local basis at each point in /C, we see that negative sectional curvature is 
ubiquitously present in our gabion state-spaces. 
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2.5.1. Recipes for constructing rules and rule fields In the context of MOR analysis, rule 
lines are easy to construct, and they have a clear algorithmic significance. 

Rules can be readily constructed via the product-sum algebraic structure of Fig. [3j by 
varying any one {'c™} while holding the others fixed. The tangents to the rule lines then 
constitute an overcomplete local basis, as depicted in Fig. |4]^C). 

More generally, rule fields can be constructed by selecting an arbitrary order (column) 
in the product-sum algebraic structure of Fig. [3j selecting arbitrary basis vectors for the 
Hilbert subspace associated with that order (equivalent to imposing an arbitrary rotation 
on the basis vectors of that column's subspace), selecting an arbitrary rank (row), selecting 
an arbitrary element of the substate of that order and rank, choosing a coordinate system 



(in the strict sense of Section 2.3) such that the selected element is one of the coordinate 
functions, and identifying a rule field on fC with the partial derivative with respect to that 
coordinate. 

From an algorithmic point of view, a rule vector is a direction in the state-space along 
which trajectories can move with great algorithmic efficiency, since only one state-space 
coordinate need be updated. 

2.5.2. The set of gabion rules is geodesically complete Whenever any two rows of the 
product-sum in Fig. [3] are algebraically degenerate, the tangent space of /C will be geometri- 
cally degenerate, yet according to the construction of the geometric rules, and in particular 
because the underlying algebraic structure is polynomial, the rules pass through curvature 
singularities without disruption of their geodesic properties, as depicted in Fig. |4]^D). 



Furthermore, it is evident that the sectional curvature (13) diverges in the neighborhood 
of a rule singularity, in consequence of the divergence of the pseudo-inverse metric g°'^ that 
appears in the projection operator Pfc as given in ([T]). Numerical experiments confirm this 
expectation. These phenomena suggest that gabion manifolds might fruitfully be analyzed 
in terms of affine algebraic varieties. The authors have not pursued this line of analysis. 

2.6. Gabion-Kahler (GK) manifolds 



We now specialize (13) to complex manifolds having a Kahlerian metric, which we will call 



gabion-Kahler (GK) manifolds. In so doing, we will adopt certain "tricks" of indexing that 



Kahlerian geometers use. We begin by writing the numerator and denominator of (13) in 
matrix notation, with all quantities still real: 

dimC dimC 

num = [(-0,,, • ^fe^) - (i/;^^, • -0 ^(/"^(f/jj • -0;,^) 

a,b,c,d=l ei/=l 
dimC 

- (tab ■ ted) + E ^tab ■ te)f^itf ' ted)] W'vW (15a) 
dimC 

den= Y [[ta-tc)^tytd)-{ta-tc)^tc-td)]^''^'^"^'' (15b) 

a,b,c,d=l 

. Now we reason as follows. S{U,V) is a real number that is independent of coordinate 
system. We are therefore free to analytically continue our coordinates, transforming (for 
example) the coordinate pair {ci,ci} via ci = xi -|- iyi and ci = xi — iyi. 

Since the sectional curvature is a geometric invariant, the (real) value of S{U, V) will not 
be altered thereby, even though the coordinates themselves are now complex. 

2.6.1. Kahlerian indexing and coordinate conventions It is evident that analytic contin- 
uation to complex coordinates treats ci and ci as independent coordinates for symbolic 
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manipulation purposes (such as partial differentiation), just as xi and yi are independent 
coordinates. It is only at the very end of a calculation, when we assign (complex) numerical 
values to ci and ci, that they can no longer be varied independently. 

It is helpful too to replace Latin indices with unbarred and barred Greek indices, in a 
convention that associates barred indices with barred coordinates (see f73^ p. 8] or [135j ) . 
Then the vector V has components {v^,v'^, . . . ,v'^ , . . .}, for example, and is represented 
in terms of partial derivatives by 



V 



dimC/2 



a=l 



dimC 



(16) 



In this convention a barred index k = dimC/2 + k, such that the index i is a shorthand for 
the integer dimC/2 + 1, and (v")* = v". 

With regard to the embedding Hilbert space Ti., we will adopt the physics convention 
that the "ket" vector ip{c) = \ip{c)) is a complex vector of dimension dim 7^/2, with the 
"bra" vector ■0(c) = (V'(c)| being the conjugate vector. 

Thus i/'(c) is a holomorphic function (also known as an "analytic function") of com- 
plex pseudocoordinates c, and ^(c) is similarly a holomorphic function of c. Defining the 
biholomorphic Kdhler potential function k(c, c) to be 



K(c,c)^i^(c).^(c) = i(^(c)|V^(c)), 
the components of the metric tensor g are given in terms of k(c, c) by 

9k(c, c) 



9al3 — 9l3a 



and gafs = g^p 



dc^dcf^ 

It immediately follows that in any holomorphic coordinate system 



0. 



9 



.al3 



5^" and g' 



.0/3 



9 



.al3 



0. 



(17) 



(18) 



(19) 



The simplest explicit example of this convention is the complex plane regarded as a 
two-dimensional Kahler manifold. Indexing its two coordinates by {1,1} yields coordinates 
{ci, cj}, which in analytic function theory are conventionally called {z, t}. The real-valued 
Kahler potential is k = ciCj/2 = zz/2., the components of the metric tensor are gah = [? o]/2 



and g"'^ = 2 [ ^ q ] , and the length element ds"^ is 

ds"^ = gab dc"" dc^ = gii dc^ dc^ + gii dc^ dc^ 

ab 



^{dz dz + dzdz) 



dzdz, (20) 



which is the usual normalization. Note the ubiquitous factors of 2 and 1/2, which require 
careful attention in practical calculations. A general scalar function on this manifold is of 
the form f{z,z), and functions of the special form /(z) are the holomorphic (or analytic) 
functions. 

In programming calculations on Kahler manifolds of larger dimension, it is helpful that 
(18) and (19) take the block-matrix forms 

and [g'^"] 



[gab] 



[ ] [9a/3] 
[9af5] [ ] 







_[g^^] [ ] . 





Q/3J 



[9a/3f 
\0] 



(21) 



where "[. . .]" is a square matrix. Environments like MATLAB and Mathematica provide 
built-in functions for block-matrix constructs of this type. We further see that in conse- 
quence of gab = gba and g^^ = (jap, which follow from (18), the individual block matrices in 
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(21 ) are Hermitian and semipositive; for this reason the submatrix [g^p] is sometimes called 
tne Hermitian metric of the complex manifold. In practical calculations it is considerably 
more efficient to work solely with the Hermitian metric and its pseudo-inverse, than with 
the larger matrix gab- 

Numerically-minded readers who are new to the literature of Kahler manifolds will 
appreciate that the above indexing conventions elegantly resolve a contradiction of our 
intuitions. On the one hand, we expect that a coordinate transformation cannot change 
the range of an index. On the other hand, we expect on physical grounds that a manifold 
described by complex coordinates will need only half the number of coordinate variables 
as the same manifold described by real coordinates. The resolution of this dilemma is in 
the block structure and symmetry properties of g, which ensure that in practical geometric 
calculations, only half the index range need be summed over. 

2.6.2. GK sectional curvature in physics bra-ket notation It is then a straightforward ex- 



ercise to write ( |15a - 15b) compactly, in the bra-ket notation of physics: 



num = ^ 



{dudutplPKldvdvip) + {dvd^i{j\Pic\duduij) - 2{dud^ij\Pic\dud^ip) 



(22a) 



den 



where the partial derivatives and the projection operator Px: = I — Px: are given in terms 
of components by 

dimC/2 

\^.^P)= Yl ^"^^IV'I^)) (^t^-) (22c) 

a=l 

dimC/2 dimC _ r.2 
a=l j3=i 

The above expressions show explicitly that bra-ket notation allows the sectional curvature 
to be computed by summing over half-ranges of coordinate indices. This notational com- 
pactness constitutes (from a Kahlerian geometry point of view) the main practical rationale 
for the bra-ket notation of the physics literature. 

2.6.3. Defining the Riemann curvature tensor Now we define the Riemann curvature ten- 
sor to be that scalar function R{A, B, C, D), defined at each point on the gabion manifold 
/C, with A, B, C, D being arbitrary vectors, such that the local sectional curvature and the 
local Riemann curvature are related by 

S{U,V) = R{U,V,U,V)/\U ^V\'' = Rabcd nVuV/|?7AF|2. (23) 

It is known (see |82| Theorem 3.8] or |135| Theorem 7.51]) that the Riemann curvature so 
defined is unique, for both real and Kahler manifolds, provided that the following index 
symmetries are imposed: 

R{A, B, C, D) = -R{B, A, C, D) = -R{A, B, D, C) = R{C, D, A, B) (24) 

which are the conventional antisymmetries of the Riemann tensor, and provided in addition 
the following identity is satisfied 

i?(A, B, C, D) + R{B, C, A, D) = R{C, A,B,D) = 0, (25) 



which is called the first Bianchi identity. 
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On real manifolds, the implicit definition (23) of the Riemann curvature in terms of the 



sectional curvature is difficult to work with, in the sense that the general expression for 
Rabcd as a function of g turns out to be too complicated to readily derive by inspection or 



manipulation of ( 13 ) . 



Fortunately, on Kahler manifolds we have the simpler definition of S{U, V) given in 



(22a 22d), from which the Riemann curvature can be read off as 



9(311,5 3 9Py,a, 



(26) 



where k{c, c) is the biholomorphic Kahler potential introduced in (17) and g^'^ is the pseudo- 
inverse of = K^^u introduced in (21). We further specify that all components of R not 
fixed by (26) vanish, save those required by the symmetries (24). Since the following 



components of R cannot be obtained from (26) by symmetry, we take them to be zero 



R. 



aficd 



R. 



aficd 



R, 



ab-yS 



R 



ab-yS 



0. 



(27) 



and we see that the resulting Kahlerian Rabcd has a block structure similar to that of the 
Kahlerian metric gab in (21 ). In consequence of this block structure, it is straightforward to 



verify that the Bianchi identity (24) is equivalent to the following Kahlerian Bianchi index 



symmetries (which we have not found explicitly given in the literature): 



R 



R. 



R 



R. 



(28) 



and which (26) respects. Thus the definition, symmetries, and identities (23-25) all are 



satisfied, and we conclude that the Kahlerian sectional curvature (22a -22d) uniquely deter- 



mines the Riemann curvature to be (26). We remark that the elegant functional simplicity 



of the Kahler-Riemann curvature tensor ( 26 ) , which in many textbooks mysteriously ap- 



pears only at the end of a long algebraic derivation, emerges quite simply and naturally in 
our Gauss-style, immersive, bra-ket derivation. 

2.6.4- The Theorema Egregium on GK manifolds Because S{U, V) and R{A, B, C, D) de- 
pend solely on the intrinsic metric g, we have thus derived — solely by analysis of sectional 
curvature — the celebrated Gauss /Riemann Theorema Egregium as it applies to Kahler man- 
ifolds. As mentioned above, we have not found in the Kahlerian geometry literature any 
similar derivation of the Riemann tensor by sectional curvature analysis. However, the 
Kahler geometry literature is so vast, that we can say with confidence only that the above 



derivation of (26) is quite different from the usual (intrinsic) derivations in the literature 
(see [135| Theorem 12.5.6] and |142| Ch. 6]). We remark also that minor misprints are 



commonplace in the mathematical literature {e.g., |120l eq. 7.11 of ch. 1] is identical to (26) 
save for an incorrect index-pair contraction). That is why we will derive some closed- form 



results in Section 2.8 against which large numerical codes can be compared. 

To better serve our MOR purposes, we have generalized R^^-yS by allowing its indices 
to range over dimC rather than dim/C, and defined the contravariant tensor g^'^ in terms 



of the matrix pseudoinverse. These extensions do not alter the functional form of (23-28), 



and they make practical numerical calculations very much easier to program, as we will see 
in Section [2?9l 

It is apparent that the generalized metric tensor g^'^ that appears in the Riemann 
curvature tensor (26) is identical to the matrix (9(Xic?k)^ that appears in the simulation 



algorithm of Step B.3 of Fig. |2] This is the first of many links that we will establish 
between Kahlerian geometry and quantum simulation. 
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2.6.5. Readings in Kdhlerian geometry Having derived the Kahlerian Riemann curvature 
tensor, we will attempt a brief survey of the Kahler geometry literature. The literature 
on Kahler geometry is comparably vast to the literature on quantum measurement and 
information, and so our review necessarily will be exceedingly sparse and subjective. 

Our development and indexing conventions in this section have paralleled the Rieman- 
nian conventions of Weinberg [195 , as extended to Kahlerian geometry by Flaherty [74J, as 
further extended to abstract notation by Martin |135j and Moroianu |142] . However, our 
analysis has been centered upon sectional curvature, rather than Riemann curvature as in 
the preceding texts. Other suitable texts include Kobayashi jl20j. Frenkel [77], Gallot 
Flanders [75], Jost |115j . Hou and Ho [lllj . and a lengthy review by Yau j200j . Many 
more textbooks and review articles on Kahler geometry exist, and it is largely a matter of 
individual taste to chose among them. 

2.7. Remarks upon holomorphic hisectional curvature 

We have seen that the sectional curvature of both real and Kahlerian gabion manifolds 



is constrained by Theorem 2.1 We now review additional (known) sectional curvature 
theorems that relate particularly to gabion-Kahler manifolds. These theorems concern a 
geometric measure called the holomorphic hisectional curvature |91j . 

We begin by remarking that physicists in particular are accustomed to thinking of if) and 
ii/j as being physically the same vector. In Kahlerian notation, the notion of "multiplying 
by i" is associated with an almost complex structure J, which is defined to be a linear map 
J, defined on every point in /C, satisfying J^V = —V for V an arbitrary tangent vector, and 
differentially smooth, that leaves the metric tensor invariant, i.e., g{U,V) = g{JU, JV). 
The effect of J upon the components {v°',v"} of an arbitrary vector V is simple: 

{v'',v"} ^ {iv'',-iv''}. (29) 



The functional form of ( 22ap2b| ) then implies the identities S{U, V) = S{JU, JV) and 
S{JU,V) = S{U,JV). 

We now consider the sum S{U,V) + S{U,JV) = S{U,V) + S{JU,V), which by defi- 
nition is called the holomorphic hisectional curvature [9.1 . The fundamental inequality of 
holomorphic hisectional curvature 

S{U,V) + S{U,JV) <0 (30) 



follows immediately from (22a 22b), because the "multiply by i" rule implies that the first 



term in the denominator cancels in the sum. We will call (30) the holomorphic hisectional 
curvature nonpositivity theorem (HBCN Theorem). 

The HBCN Theorem is a long-known result 1*911 that applies in general to any Kahler 
manifold that is a complex submanifold of a Euclidean space. Physically speaking, if a 
given Kahlerian section {[/, V} has positive curvature, then the "rotated" section {[/, JV} 
will have negative curvature. It is readily shown [21] that the HBCN Theorem implies the 
nonpositivity of the eigenvalues of the Ricci tensor, and therefore, the nonpositivity of the 
scalar curvature. 

As a technical point, the HBCN Theorem applies only to Kahler manifolds that have 
(or can be given) a complex embedding in a larger Euclidean manifold, which is the case 
of greatest interest in quantum MOR applications. The HBCN Theorem does not apply 
to Kahler manifolds that have no complex Euclidean embedding. For example, Kahler 
manifolds having a Fubini-Study metric can (and sometimes do) exhibit positive scalar 
curvature. This is incompatible with the HBCN Theorem, and these manifolds therefore 
have no complex Euclidean embedding. 
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The functional form of (22a) allows us to immediately extend the HBCN Theorem to 



encompass the directed sectional curvature S{U, V, h) as follows: 

Theorem 2.2 (directed extension of the HBCN Theorem). At an arbitrary point on a 
Kdhler manifold fC that has a complex embedding within a Euclidean space TC, let U and V 
be arbitrary tangent vectors and let h be an arbitrary unit vector normal to /C. Then the 
directed sectional curvature satisfies S{U, V, fi) + S(U, JV, h) < 0. 



2.7.1. Relation of Theorem 2.1 to the HBCN Theorem 



from both the HBCN Theorem and its directed extension Theorem 2.2 



Our Theorem 2.1 is a different result 
Most obviously, 



Theorem 2.1 applies to all manifolds that possess one or more rule fields, whether they 



are Kahlerian or not, while both the HBCN Theorem and our Theorem 2.2 apply solely to 
Kahler manifolds. 



Even on Kahler manifolds, the HBCN Theorem and Theorem 2.1 have substantially 



differing implications. It follows from the defining equation of a rule field (14) that if W is 



a rule field, then so is JW, hence if 1/ is a rule vector, then so is JV . It then follows from 



Theorem 2.1 that S{U.,V) and S{U,JV) are both nonpositive, which is a strictly stronger 
condition than the nonpositivity of their sum that is implied by the HBCN Theorem. 



It further follows immediately from (15a -15b) that if V is a rule vector, then S{U,V) = 
S{U.,JV). More generally, if U and V are both rule vectors, then S{U.,V) = S{U.,JV) = 
S{JU,V) = S{JU, JV), and all of these sectional curvatures are nonpositive. From this 
"gabionic" point of view, we see that on gabion-Kahler manifolds the rule vectors U and JU 
are effectively the same vector, i.e. rule vector phases are irrelevant to sectional curvature 
properties, as is natural in quantum mechanical analysis. 

The ubiquity of nonpositive sectional curvature on GK manifolds therefore can be viewed 
as arising from the confluence of rule structure and complex structure, both of which are 
associated with strong theorems that imply negative sectional curvature. 

2.7.2. Practical implications of sectional curvature theorems From an MOR point of view. 



the nonpositive sectional curvature implied by Theorem 2.1 can be regarded as helping 



to ensure the robustness of MOR on both real and complex manifolds, while the HBCN 
Theorem helps makes MOR even more robust on complex manifolds. This is one of two 
fundamental reasons why quantum MOR can be regarded as intrinsically easier than classical 
MOR (the other reason being the algorithmic resources that are provided by the Theorema 
Egregium, as we will discuss in the following section). 

We remark that not all sectional curvatures S{U,V) are negative on GK manifolds. 
In small-dimension GK manifolds we have constructed analytical examples of positive- 
curvature sections — with some trouble because neither U nor V can be rule vectors — and in 
large-dimension gabions numerical searches find them. According to our present (limited) 
understanding, these positive-curvature Kahlerian sections seem rather artificial, and so we 
will not discuss them further. It is possible that their significance has eluded us. 

2.8. Analytic gold standards for GK curvature calculations 

The Riemann tensor specifies the Gaussian curvature of every section of /C, and on large- 
dimension MOR manifolds the Riemann tensor therefore carries a vast amount of informa- 
tion. To compress this information (among other purposes) it is conventional to condense 
the four-index Riemann tensor Rabcd into the two-index Ricci tensor Rab by 



dimC 



(31) 



c,d=l 
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which in Kahlerian index notation can be written 

dimC dimC/2 dimC 

c,d=l 7=1 S=l 

Further compression is achieved by the scalar curvature R defined by 

dimC/2 dimC 



(32) 



dimC 

a, 6=1 



2E E^"%- 



(33) 



a=l 



13=1 



Caveat: various authors entertain diverse conventions regarding which indices should be 
contracted to obtain the Ricci tensor, and the factors of -1 and 2 that appear above are 
commonly associated with minor errors and imprecisions. 

There is one sectional curvature convention, however, that is universal: the directed 
sectional curvature on a unit hypersphere 5" (the ordinary sphere embedded in R^ being 
S^) satisfies S{U,V,h) = 1 for all linearly independent U and V. In a locally orthonormal 
basis there are n(n — 1) such pairs; it follows that on a unit hypersphere the eigenvalues 
of the Ricci tensor are all n — 1, and the scalar curvature itself is therefore R = n{n — 1). 
This result is a useful "gold standard" for testing symbolic and numerical calculations on 
real manifolds. 

To construct a similar gold standard for testing symbolic and numerical curvature cal- 
culations on Kahlerian manifolds, it is convenient to consider the manifold of rank-one, 
order-n gabion states (see Fig.|3]). We allow the spin quantum numbers {ji : i E l,n} asso- 
ciated with successive product subspaces to vary independently. Then by a straightforward 
(but not short) calculation, it can be shown that the scalar curvature of a general rank one, 
order-n product state is 



R 




k,m=l 



for m ^ k 
otherwise 



(34) 



where k is the Kahler potential function of (17). See the following section for a summary 



of the assorted algebraic techniques used to obtain this result. 

To the best of our knowledge, this general analytic result has not previously appeared 
in the literature. For the simple example manifold of (|3|, we have order n = 2 and ji = 
j2 = 1/2, so the above yields R = —4/k = —8/{ip\'ilj) = —8/{il}-ip), in agreement with the 



automated "push-button" analysis of Section |2.3 

2.9. The Riemann-Kahler curvature of Slater determinants 

We now have all the tools we need to compute the scalar Riemann curvature of Slater 
determinant states, which are the main state-space of quantum chemistry |54| \80\ I184j . The 
strategy of the calculation is straightforward, and although the details are lengthy, the final 
result is simple. 

We begin by considering, as the simplest example having nontrivial curvature, the quan- 
tum states that are obtained by antisymmetrizing the outer products of a rank-one order-two 
product-sum state of spin j = 3/2 (see Fig. Isl) 
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In the language of quantum chemistry, we can equivalently regard i/'(c) as the variational 
state-space of the (unnormalized) Slater determinant states of two electrons "a" and "6" 
occupying linear combinations of four orbitals. 

With i/'(c) given, we can compute its scalar curvature R by either numerical or ana- 
lytic means. To compute R numerically, we can simply set the eight pseudo-coordinates 
c = {c^,c^; i = 1,4} to any desired value (thereby choosing the point in state-space at 



which R is to be evaluated) then compute first the Kahler potential k from (17), then the 



Riemann curvature R from (26), and (31-33), evaluating the pseudo-inverse metric g^^'^ of 



(26) numerically. 



Empirically we find that for our simple example (35) these numerical calculations in- 
variably yield a Riemann scalar curvature of i? = — 8/k (to machine precision) for all values 
of the pseudo-coordinates {Cq,c^. 

The simplicity of the numerical result motivates and guides the following analytic eval- 
uation of R in closed form. We designate by ■j/'q the point in the state-space ■?/'(c) at which 
the Riemann curvature is to be evaluated. We write ipQ in the following standard form by 
an appropriate choice of basis vectors 
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An arbitrary state in a neighborhood of Vo can be written as an explicit holomorphic 
function of five complex coordinates {c^,c^,c^,cl,c'^} as follows: 
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That the above variety is in fact equivalent to (35) can be demonstrated by the "push- 
button" method of computing their respective Groebner bases |53.j and verifying that these 
bases generate the same homogeneous ideal in the variables { 

In the transition from (35) to (37), four pseudo-coordinates {Ca,c^,cl,c^} have disap- 
peared and been replaced by the "0" and "1" entries. These are physically accounted as 
follows: and have been merged into the overall (complex) normalization by a simple 
rescaling; and cf have tangent vectors that vanish at tpo, such that they do not induce 
coordinate charts on i/'(c) at -00 and can safely be dropped. In numerical calculations, their 
tangent vectors at tpo are in the null-space of g^^ and 5"^. 

The remainder of the analytic calculation is straightforward. The metric tensor compo- 
nents are given from (37) by (17 18), and when evaluated at i/jg, yield a 5x5 matrix 



that is diagonal and nonsingular. Computing the inverse g'^'^ is therefore trivial. The Rie- 
mann components RajsyS given by (26) and the Ricci components R^p are given by (32 ). 



^Specifically, the Grobner basis of both (351 and (37l is found to be {t/'i, i/Je, i/'ii, ^Aie, ''/'2 + tp5,ip3 + 
^9,''/'7 + i/'io,''/'i2 + V'i5,''/'4 + V'i3,''/'8 + V'14, V'lo'/'is + ''/'5V'i5 — i'-ii'ii] ■ The disadvantage of the Grobner basis 
is that it obscures the symmetries of the parametric representation ( |35[ ); for example the rule structure of 
the manifold is not evident. Furthermore, implicit representations like this one are poorly suited to sectional 
and Riemannian curvature calculations. 
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Finally, the scalar Riemann curvature given by (33) is found to be i? = — 8/k, in accord 
with the numerical result. 

This reasoning is readily generalized. Upon antisymmetrizing a general rank-one, order- 
n product-sum state (see Fig. |3]) under the exchange of all pairs of spins — thus converting it 
to a single n-particle Slater determinant — and evaluating the metric and Riemann tensors 



in the diagonal basis of (36-37), the dimensionality of the Kahler manifold is evidently 



dim/C = 2(1 -|- n(norb — n)), and the scalar curvature is found in closed analytic form to be 



^ _ _2 ^ fn(n - l)(norb - n)(norb - n - 1) for > n ^^^^ 
n 1 undefined because ip{c) = for riorb < n • 



Here norb = 2 j -|- 1 is the dimension of the individual states in the Slater product-sum, with 
the mnemonic "orb" referring to "orbitals" in token of the practical use of these states in 
quantum chemistry. Since the scalar curvature is a geometric invariant, this result holds in 
any basis. 

The curvature of a Slater determinant manifold having a Fubini-Study metric, i.e., 
having a Kahler potential k = 1/2 log('0|i/') instead of k = l/2{ip\ip), can similarly be 



calculated from (37) and its multidimensional generalizations. The result is 



^ _ 1 4nnorb(raorb - n) for norb > n ^^^^ 

[undefined because ■j/'(c) = for norb < ■ 

Physically speaking, the Fubini-Study metric describes a manifold of normalized states in 
which a phase rotation of 6(f) radians has a path length of zero rather than 6(j). Furthermore, 
it is a straightforward (but lengthy) algebraic exercise to verify during the course of the 
above calculation that g^j^ cx R^p: Slater determinant manifolds having a Fubini-Study 
metric are Ricci- Einstein manifolds [13] • An immediate consequence is that Slater deter- 
minant manifolds are solitons under Ricci flow [39]. We note that Ricci soliton manifolds 
are of central importance to the mathematical community, because they represent (in sense 
that can be made precise) the unique "smoothest" manifolds of a given topological class. In 
recent years this idea has been central to the Ricci flow [50] proofs of several long-standing 
topological conjectures. 

To the best of our knowledge, the above algebraic/geometric properties of Slater de- 
terminants have never been noted in the literature. This gap is noteworthy, in view of 
the central role that Slater determinants play in chemistry and condensed-matter physics, 
and the similarly central role of Kahlerian algebraic geometry in numerous branches of 
mathematics. It reflects what the National Research Council has called 11491 



"[...the anomaly that] although theoretical chemists understand sophisticated 
mathematics and make heavy use of the mathematical literature, they have typ- 
ically not involved mathematicians directly in either the development of models 
or algorithms or the derivation of formal properties of equations and solutions. 
In fact, theoretical chemists have become accustomed to self-reliance in mathe- 
matics." 

A central objective of this article's geometric approach to quantum simulation is to more 
closely link the formal mathematical tools of algebraic geometry to practical problems in 
quantum simulation, and thereby, to help ensure that the emerging discipline of quantum 
system engineering is not needlessly "self-reliant in mathematics." 
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2.10. Slater determinants are Grassmannian GK (GGK) manifolds 

We thank our colleague Joshua Kantor for directing our attention to the facts that the 
Slater determinants associated with n particles distributed among norb orbitals have a 
natural isomorphism to the Grassmannian manifolds classified as G{n, Uorh) |100| . that 
Grassmannian manifolds have a known presentation as projective algebraic varieties via the 
Pliicker embedding [7T] , and that the Pliicker embedding is possessed of isometries having 
a Lie group structure that is known to be compatible with a Kahler-Einstein metric [13] • 
From this algebraic geometry point of view, the example of a Slater determinant that 



we present in (35), having n = 2 and norb = 4, can be identified with what Harris' classic 



textbook Algebraic Geometry calls MOO- p. 65] "the first nontrivial Grassmannian — the first 
that is not a projective space — [namely] G(2,4)." Further investigation of this convergence 
of algebraic geometry with quantum chemistry is in progress. 

2.11. Practical curvature calculations for QMOR on GK manifolds 

Now our appetite is whetted to ask even more difficult questions: what are the curvature 
properties of higher-rank GK manifolds? And what are the implications of these curvature 
properties for quantum simulation? We will turn to numerical experiments to learn more 
about these issues. 

As we begin numerical curvature calculations on Kahler state-spaces of increasingly 
large dimension, we first consider what we may expect to learn from these calculations that 
would have practical implications for MOR. 

It is known [78i ,144j that at each point on IC we can construct locally Euclidean Riemann 
normal coordinates such that the local curvature tensor has the expansion 

dim?C 

gab = Sab-Y. I RiO)cadbyV (40) 
c,d=l 

and the local volume element dV oc | det g |^/^ therefore has the expansion 

(dim/C \ 
1- lR{0)aby''y'] . (41) 
a,b=l I 

We see that the square roots of the eigenvalues of the Ricci tensor determine the length scale 
over which curvature effects dominate the volume of the state-space. Physically speaking, 
this sets the length scale over which gabion petals are locally Euclidean. We further see that 
negative curvature is associated with an exponential "flowering" of the state-space volume, 
in accord with the beautiful knitted representations of hyperbolic spaces by Taimina and 
Henderson |1L)9] . 

We note that since the Ricci tensor eigenvalues are geometric invariants (meaning specif- 
ically, the eigenvalues of the mixed Ricci tensor R°^ = '^^g"'^Rcb are geometric invariants), 
we can compute them in any coordinate system we please. We further note that our adoption 
of a over-complete Kahlerian basis creates no anomalies in the Ricci eigenvalue distribution, 
because the extra eigenvalues of R"^ = g"''^ R^/s vanish identically, due to the projective 



definition (21) oi g°''^. 

Our first goal, therefore, will be to compute the Ricci tensor eigenvalues for high-rank 
GK manifolds, expecting thereby to gain quantitative insight into the "flowering" gabion 
geometry that is depicted in Fig. |4]^H). 

As our first test case, we will consider the rank-6, order-9, spin-1/2 GK manifold. This 
manifold has (real) dimension dim/C = 2x9x(6-|-l) = 126 and it is embedded in a 
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Figure 5: Typical Ricci tensor eigenvalues for gabion-Kahler manifolds. 

Points on the gabion-Kahler manifold are randomly selected by first 
randomly generating an independent normalized product state for each 
gabion rank {i.e., each row of Fig.|3j), then summing the ranks, then nor- 
malizing the final state. The sparsely dotted line contains typical Ricci 
eigenvalues of the spin-1/2, rank-9, order-6 gabion, the densely dotted 
line contains the eigenvalues of the spin-1/2, rank-9, order-10 gabion. 
The dashed line is an empirical "rank x order" estimate of the mean Ricci 
eigenvalue of gabion-Kahler manifolds having large-co dimension. 



Hilbert space of (real) dimension dimTi = 2x2^ = 128. The gabion state-space therefore 
"almost" fills the Hilbert state-space, since only codim /C = dimTC — dim/C = 2 dimensions 
are missing. 

Choosing a random point in /C and computing the Ricci eigenvalues numerically yields 
the typical eigenvalues shown in Fig. [5} We take the square root of the largest eigenvalue 
to be (roughly) the linear extent of a gabion petal; these petals evidently have an extent of 
order 0.01. 

We emphasize that IC is large enough to contain exponentially many such petals. See 
Nielsen and Chuang for a quantitative analysis |146| Section 4.54], noting that the "patches" 
of Nielsen and Chuang are broadly equivalent to our petals and therefore, their Fig. 4.18 
is broadly equivalent to our Fig. |4]|^H). This mathematically justifies our physical picture 
of a gabion as a geometric "flower" having exponentially many petals. Later on, in Sec- 



tions 4.6.4-4.6.6 we will rigorously count the number of petals using ideas from coding 
theory. 

2.12. Numerical results for projective QMOR onto GK manifolds 

Suppose we generate a random (normalized) quantum ■j/jq, and numerically search for a 
highest-fidelity projection of i/^q onto IC. We will call this high-fidelity image point i/)^. 
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Based on our geometric analysis so far, we have two strong expectations relating to quan- 
tum MOR by projection onto /C. We expect first, that projections exist for which — "001 
is no greater than ~ 0.01, since at greater separations the negative curvature of /C is great 
enough to ensure a better solution, via the mechanism of Fig.[4j^G). Second, we expect that 
the numerical search for tpj^^ will be well-conditioned. That is, it will robustly converge to 
a high-fidelity representation, despite the exponentially convoluted geometry of /C, with- 
out becoming trapped in local minima (because every local section that contains a rule is 
negatively curved). 

In numerical trials both expectations are fulfilled. The achieved fidelity is in good accord 
with the geometric expectation of ~0.01: the median value of {if^ji^ — in a trial of 100 
projective reductions was 0.005, and the maximal value was 0.025. 

We computed the reduced-order GK representation tpj^ by a simple gradient search. 
Specifically, we integrated to convergence the dynamical equation of Step B.3 of Fig. |2] 
with the potential cp replaced by 

'^-^^(c),a-(^o-0^(c)). (42) 

Upon conversion to Kahler index notation, the resulting reduction equation is manifestly 
geometrically covari ant : 



dimC r ~i / \ 

0(c) 



/3=i 



(^o-^(c)) 



(43) 



From an algorithmic point of view, this reduction equation can be regarded as a downhill 
gradient search for a reduced-order representation tp^^ = lim^^oo ■'/'(c(t)) of ipQ. The search 
evidently converges when ipK^ is directly "underneath" ■j/>q, i.e, when Pfc{'ipi^)(^ipQ — 'ipi^^ = 0, 
where we recall that Pic{ipfc) projects vectors in 7i onto the tangent space of /C at ipj^. As an 
alternative to gradient search, we note that Beylkin and Mohlenkamp sec. 3.1] describe 
an alternating least-squares algorithm that in their hands gives excellent results, but we 



have not ourselves tried this method. In Section 4.6.1 we derive the above equation (43) 
using the language (and nomenclature) of compressive sampling (CS), and we describe the 
numerical calculations in greater detail. 

No attempt was made to optimize the efficiency of the gradient search; instead we simply 
reused the existing dynamical code implementing Step B.3 of of Fig. [2] (this code exploits the 
gabion algebraic structure of ip{c) to evaluate ( |43[ ) efficiently). This duplication of internal 
algorithmic structure again illustrates the intimate relation between trajectory calculations 
and geometric calculations. 

No gross failures of convergence, such as might be expected from trapping of the tra- 



jectory induced by (43) on a distant gabion "petal," were observed in this (or any) of 
our numerical trials. Our geometric analysis explains this robustness as originating in the 
negative sectional curvature of the ruled net of the gabion state-space /C, as depicted in 
Fig.gG). 

We now increase the order of the gabion to 10, keeping the rank at 9. This increases the 
dimensionality of the Hilbert space to 2 x 2^*^ = 2048, and the dimensionality of the Kahlerian 
gabion to dim/C = 2 x 9 x (10 1) = 188. With codim/C = 1048 - 188 = 1860, our model 
order reduction is now discarding ~ 90% of the dimensions of the larger Hilbert space. Thus, 
in deliberate contrast to the previous example, the MOR is now fairly aggressive. Typical 
resulting Ricci eigenvalues are shown in Fig. [Sj 
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A pronounced flattening of the eigenvalue distribution is evident in this large-codim 
example. The plotted straight line is simply rank x order, which seems empirically to 
describe the average Ricci eigenvalue in this particular case, and also in many other trials 



that we have run. This rule-of-thumb is simply the analytic result for rank-1 curvature (34) 
multiplied by the rank. 

There is at present no analytic theory that justifies this rule-of-thumb, or explains the 
observed flattening of the eigenvalue distribution. Obviously, such a theory would be wel- 
come. We remark that to the best of our knowledge, the Ricci eigenvalues of Fig. [5] having 
dim/C = 188, are the largest-dimension curvature eigenvalues ever numerically computed. 



We then calculated projective MOR approximations by integrating (43) as before. We 
generated our targets tpQ by randomly selecting target points on /C, and moving a distance 
of 0.05 along a random vector h perpendicular to IC. 

Again, robust convergence to high-fidelity reduced-order representations was observed. 
Of 100 trials, in 98 cases the separation distance was precisely 0.050, representing conver- 
gence to the correct "petal." In the remaining two cases the separation distances were 
0.123 and 0.120, representing sporadic convergence to a wrong-but-nearby "petal." Thus 
even wrong-petal convergence yielded a high-fidelity representation. This robustness can 
again be ascribed to the negative sectional curvature of the state-space manifold's ruled 
net, according to the geometric mechanism depicted in Fig. QG). 

2.13. Avenues for research in geometric quantum mechanics 

To begin, the results of the preceding section are understood only qualitatively, and rigorous 
bounds on projective fidelity and robustness would be very welcome. How much of the 
preceding section can be understood solely as a consequence of the known sectional curvature 
properties of the GK manifolds? 

Although our presentation focuses on practical applications of quantum MOR, a broad 
class of fundamental physics questions can be given a geometric interpretation. We tem- 
porarily adopt the point of view of geometric quantum mechanics in which — as reviewed in 
Section |1.3| — the manifold /C is regarded as the "real" arena on which physics takes place. 

We begin by noting that that the rule-field equation ( |14[ ) specifies the dimensionality of 
Hilbert space to be the number of (linearly independent) scalar rule-fields that the (postu- 
lated) "real" GK manifold of geometric quantum mechanics supports. The question then 
arises, what determines the number of these embedding state-space fields? The present 
article suggests that this question is best investigated from a blended informatic-algebraic- 
geometric point of view. 

As another example, suppose p and q are arbitrary linear Hermitian operators in the 
embedding Hilbert space 7i. At a given point of the Kahler manifold /C, specified by a state 
\jp{c)), we can construct tangent vectors Vq and Vp whose components are 

dimC o 

< = E ^"'^ ^(^(c)klV'(c)) and = {v^r (44) 

0=1 ^ 



and similarly for Vp. With this normalization we have g{Vq,Vq) = {iplqP/cqltp) . Physically 
speaking, Vq defines a Kahlerian velocity field along which the projected dynamical equation 
d\ijj{t))/dt = Picqlipit)) moves state trajectories. 

The sectional curvature S(yp,Vq) then can be regarded as a fundamental property of 
the "true" physical manifold IC. The Theorema Egregium guarantees that the sectional 
curvature is an intrinsic property of /C, and our physical intuition suggests that it should 
therefore be measurable. 
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A host of fundamental questions then arise quite naturally, that intimately unite physics 
and mathematics in the context of quantum simulation. If the quantum sectional curva- 
ture is physically measurable — whether in reality or within projective simulations — then 
by what kinds of experiment? Are these experiments practical in our real world? Is there 
any experimental evidence already at-hand to indicate that quantum sectional curvature 
vanishes in the physical world? If so, to what precision has this been verified? 

It is apparent that quantum chemists can apply a reverse strategy to create a geomet- 
ric context for assessing the fidelity of chemical quantum simulations. Set p to be the 
kinetic energy of the electrons of a molecule, and q to be the potential energy (including 
perturbations due to applied potentials), and set the gabion state-space to be a sum of 
Slater determinants. Then the vanishing of the sectional curvature S(yp, Vq) indicates that 
the state-space section generated by {p, q} is Euclidean, as is (presumably?) desirable in 
chemical simulations, most particularly, density functional theory (DFT) simulations. 



And finally, to anticipate, in Section 3.2.1 we will discuss how the Theorema Dilectum 



manifests itself upon GK state-spaces. This will turn out to raise thorny issues of how 
causality works in geometric quantum mechanics. For the present we will say no more about 
these difficult fundamental questions, instead referring the reader to Leggett's recent review 
|126] . and we return instead to our central topic of practical quantum spin simulations. 

2.14- Summary of the geometric analysis 

We have seen that the ruled net of a GK manifold, with its associated nonpositive sectional 
curvature, plays a geometric role in quantum simulation that is broadly similar to the role 
of polytope convexity in linear programming, namely, the role of providing geometric foun- 
dations for developing efficient and robust algorithms. However, understanding polytope 
convexity has been a slow process, and this research is still being pursued along multiple 
mathematical lines that include algebra, geometry, information theory, and their numerous 
hybridizations. We foresee that the role of algebraic geometry in classical and quantum 
simulations will be elucidated along similarly multidisciplinary lines, and we are conscious 
that the analysis presented here has barely begun this enterprise. 

Even in its present early stage, geometric analysis provides grounds for confidence that 
quantum model order reduction and trajectory simulation can be achieved with high ef- 
ficiency, fidelity, and robustness, provided that some physical mechanism can be found to 
compress quantum trajectories onto the petals of gabion-type state-spaces. That needed 
physical mechanism is, of course, the Theorema Dilectum, which will be the topic of the 
section that follows. 

3. Designing and Implementing Large-Scale Quantum Simulations 

Our goal in this section is to join the geometric ideas and theorems of the preceding section 
with the well-known ideas and established design principles of linear quantum mechanics. 

We have reason to worry that perhaps very few principles of linear quantum mechanics 
will survive the transition to reduced-order quantum MOR mechanics. After all, we have 
seen that gabion-Kahler (GK) state-spaces are strongly curved, so that when we project 
quantum trajectories onto them, we (seemingly) dispense with all of the mathematical 
properties of Hilbert space that depend upon its linearity. Furthermore, depending upon 
the degree of the model order reduction imposed, the GK projection of QMOR mechanics 
may even discard all but an exponentially small fraction of the dimensions of the embedding 
Hilbert space. 

We seek, therefore, to establish that the principles of linear quantum mechanics hold 
true in QMOR simulations, and in particular, to establish precisely the mechanisms by 
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which they can hold true. 

3.1. Organization and nomenclature of the presentation 

For convenience, whenever we derive a result that is novel or expressed in a new form, 
we state it as a formal (numbered) design rule that is accompanied (as needed) by formal 
(numbered) definitions. On the other hand, whenever a result is unsurprising, or can 
be found in the literature, and is obtained from previously given equations by standard 
manipulations and reductions, we outline the derivation but omit details. 
We define the subject of our analysis to be QMOR mechanics: 

Definition 3.1. QMOR mechanics is the mechanics of a physical system simulated accord- 
ing to the orthodox principles of linear quantum mechanics, as modified by projection onto a 
lower- dimension manifold having a Kdhler geometry, for purposes of quantum model order 
reduction (QMOR), as concretely embodied in the algorithms and algebraic structures of 
Figs, m 



We seek to construct recipes by which QMOR mechanics simulates linear quantum mechan- 



ics as closely as feasible. To recapitulate the objectives of Section 1.1.1 our analysis seeks 
to be orthodox in its respect for linear quantum mechanics, operational in the traceability of 
its predictions to measurement processes, and reductive in the sense that its principles are 
summarized by closed-form analytic design rules. Our analysis strives also to be synoptic 
in the sense that whenever we choose between equivalent analysis formalisms, we state a 
rationale for our choice, and we note the practical consequences of alternative choices. 

3.2. QMOR respects the principles of quantum mechanics 

We begin by establishing that the mathematical structure of QMOR mechanics is sufficiently 
rich to respect the following principles of quantum and classical mechanics: 



the causal invariance of the Theorema Dilectum is respected (Sec. 3.2. 1[ ) 



the entropy of systems in thermodynamic equilibrium is respected (Sec. 3.2.2) 



• the principles of classical linear control theory are respected (Sec. 3.2.3). 

• the quantum limits to measurement noise and back-action are respected (Sec. 3.2.6). 

If QMOR mechanics did not respect these principles, it would scarcely be useful for simulat- 
ing real-world quantum systems. Conversely, by stating these principles in a quantitative, 
we gain a fairly clear picture of the path that our QMOR analysis needs to take. 

3.2.1. QMOR respects the Theorema Dilectum We now state the Theorema Dilectum in 
an algebraic form that is well-adapted to the formal quantum MOR algorithm of Fig. [T] 
Our mathematical analysis parallels that of Nielsen and Chuang [146i see their Theorem 8.2 
in Section 8.2], whose analysis in turn derives from a 1975 theorem of Choi 08] (see Sec- 
tion 



1.3.4). 



We suppose that at the end of step n the QMOR simulation algorithm of Fig. [T] has 
computed a wave-function \ipn)- For simplicity, we further suppose that precisely one mea- 
surement operator pair {M(_,_), M(_)} acts during the subsequent time-step. These operators 
satisfy the normalization condition 

M(+)M/^) + M(_)M/_, =/ (45) 

where / is the identity operator. In the absence of GK projection the post-timestep density 
matrix pn+i is readily shown to be 

Pn+l = M^+)PnMl^ + M(_)P„M(t) = C[pn] (46) 
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where pn = \tpn){'^n\- This expression imphcitly defines the well-known linear superoperator 
C to be that linear operation on (Hermitian) matrices that takes pn Pn+i- The existence 
and strict positivity of this linear map is one of the main defining characteristics of linear 



quantum mechanics (as reviewed in Section 1.3.4). 

Now we ask "What mathematical operations upon {Af( 
Theorema Dilectum in its algebraic form due to Choi [48 
operations are of the general form 



M, 



(+) 



U 



M, 



(+) 
(-) 



} leave C invariant?" The 
states that all such invariance 



(47) 



where U is an arbitrary 2x2 unitary matrix of complex numbers (i.e. a matrix acting on 
the linear space of measurement operators, not the Hilbert space of \ipn), such that the 
matrix elements of U are c- numbers). This is the sole general mathematical invariance of 
the first two steps of the simulation algorithm of Fig. [T| and so (from the QMOR point of 
view) it is the most fundamental invariance of linear quantum mechanics. 



In Section 3.3.6 , we will establish that [/-transform invariance enforces physical causality. 
We are thereby motivated to ask "How is the [/-transform invariance of the Theorema 
Dilectum affected by GK projection?" According to the algorithm of Fig. [T] GK projection 



modifies (46 1 to 



Pn+l = (PK)nM(+)|V'n)(^„|Mi)(PK)r 



+ iP^)nM^_^\iJn){i^n\Ml^{P^)n 



K:;nM(+)|V'n) 



(48) 



where we recall that {PK)n = -f^d^n)) was defined in ([7| as the operator that projects onto 
the local tangent space of the GK manifold at iV'n)- It is easy to check that tr[/3„-|_i] = 1 



{i.e., probability is conserved), and that above expression reduces to the linear result (46) 



whenever the commutators [{Px)n,M^ 



(+). 



and [(P;c)n,M(_)] 



both vanish. 



As it stands, the nonlinear projective evolution (48) is not [/-transform invariant, and 



hence it does not always respect the Theorema Dilectum. We therefore introduce the further 
assumption (which our QMOR simulations will always respect) that both M(+) and M(_) are 
near to being multiples of the identity operator. We quantify "near" by introducing an 
artificial expansion parameter e such that 



M(+) =cl + e dM^+^ 5M(_) = si + e 5M(_) 



(49) 



where 



c = tr M(+) / dim M(+) s = tr M(+) / dim M, 



{+) 



(50) 



Equations (49-50) uniquely define the products e(^M(+) and e5M(_y Therefore we can 



uniquely define e, and thereby uniquely define (5M(+) and 5M^_) too, as follows: 



tr[(e (5M(+))t (e <5M(+)) + (e (5M(_))t (e (^M(_))]/ dim M(, 



+) 



Equations (45-51) then imply the exact normalization relations 

|2 

S| -|- |fc 

tr[(5M(+))t (5M(+,) + (5M,_))^ 



IcP + IsP + lei 



1 



dimM, 



(+) 



(51) 

(52) 
(53) 
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,(+) and (^M(_) to be 0(1). 



and so we can regard all operator products involving 5M(. 

In aggregate, our definitions and normalizations ensure that Pn+i as given by ( |48| ) has 
a well-defined power series expansion in e; it is therefore a straightforward (but not short) 
calculation to verify that this expansion can be written in the form: 



Pn+l = {PK)nC[\i>n){^n\] {P^)r 
+ \i'n){i^n\ tr [(Pk)„£[ 



3 1 



(54) 



where Pic = I — P^.- We note that the leading terms are determined solely by C and 
and therefore are invariant under the JJ-transform (47) of the Theorema Dilectum. 



Furthermore, it can be shown that the value of the small parameter e is itself invariant 
under the [/-transform, and so is the sum |cp -|- And finally, we have calculated the 
(exact) c- and s-dependence of the O(e^) terms. 



We will establish in Section [3. 3| that in the continuum limit of infinitesimally small time 
step intervals 5t, physical quantities (for example, relaxation rates) are 0{e'^/6t). Physically 
this means that the 0{e'^) terms in (54) are negligible in the continuum limit, provided the 



technical conditions |c| > and |s| > are satisfied (these technical conditions provide the 
rationale for calculating the c- and s-dependence of the 0{e^) terms in (54) ). These results 
motivate us to adopt from Carlton Caves |36] the following definition: 



Definition 3.2. Measurement operations satisfying e <C |c| and e <C |s| are called measure- 
ment operations of the first class (or sometimes first-class measurements J. 



The result (54) then expresses the following fundamental design rule of QMOR mechanics: 



Design Rule 3.1. In the continuum limit, quantum trajectory simulations of first-class 
measurement processes, as projected onto state-spaces of gabion-Kdhler type ( QMOR simu- 
lations for short), respect the unitary invariance of the Theorema Dilectum. 

In physical terms, first-class measurements are characterized in the continuum limit by 
stochastic drift and diffusion processes on GK manifolds, rather than quantum jumps. To 
ensure that the Theorema Dilectum is respected, QMOR mechanics must therefore simulate 
all quantum systems wholly in terms of drift and diffusion processes upon strongly-curved, 
low-dimension state-space manifolds of gabion-Kahler type. 

This low-dimension curved-geometry stochastic description of QMOR mechanics thus 
contrasts sharply with the high-dimension linear-geometry "jump-oriented" description of 
quantum mechanics that is commonly given in textbooks; that is why many of our derivation 
methods and design rules are novel. 

Turning our attention briefly to the broader context of geometric quantum mechanics 



(see Sections 1.3.5 and 2.13), we emphasize that we do not regard issues of causality in 
geometric quantum mechanics as settled by the projective generalization of the Theorema 



Dilectum given in (54). We will not discuss this topic further, partly for reasons of space. 



but mainly because we regard it as being an exceptionally difficult and subtle topic that is 
intimately bound-up with the question of the existence (or nonexistence) of quantum field 
theory in geometric quantum mechanics. 

3.2.2. QMOR respects thermal equilibrium Our recipes for simulating contact with thermal 
reservoirs in QMOR mechanics yield an algebraic result (not previously known) that holds 
exactly even in linear quantum mechanics, and that can be verified without reference to 
drift and diffusion equations. We state and prove this result now, so that it can provide a 
well-defined mathematical target for our subsequent QMOR analysis. 
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We begin with a brief summary of coherent states, referring the reader to classic text- 
books, such as those by Gardiner ^85j, Gottfried [93j, Klauder and Skagerstam |118] . Perelo- 
mov [153] ■ Rose |169j . and Wigner |198] . for details. We will mainly follow Gottfried's 
notation. 

We start by identifying a spin=j state having z-axis quantum number m with the ket- 
vector |j, m). Then a coherent state \x) associated with a unit- vector spin direction x is by 
definition \x) = D((j), 6, 0)|j, j)) where the rotation operator D that carries t = (0, 0, 1) into 
X = (sin cos (j), sin 9 sin </>, cos 6) is 

D{(t), e, tP) = e-*<^«3e-*^"2g-i^.3 (55) 



The rotation operators are well understood. In particular, an identity due to Wigner |93| 
1169) gives {j,m\x) in closed form as 

{j,m,\x) = Dl.{^,e,0)=( y^V^-<^(cosie)^+'"(sini^)^-'" (56) 

It follows that {x\x) = 1 and (a;|s|a;) = jx. It is well known [851 1153| I154| I162j (and not 
hard to show from (56) ) that a resolution of the identity operator I is 

, 2i + 1 



/ d^x\x){x\. (57) 

Jatt 



The Q -representation and P -representation of a Hermitian operator p are then defined 
(following Perelomov's conventions [153] ) as 

Q{x\p) = {x\p\x), (58) 



P=^4^ / d^'x P{x\p)\x){x\. (59) 



4tt 



Given an arbitrary Hermitian operator p, it is known that in general a P-representation 
P{x\p) can always be constructed. In brief, the construction is as follows: from the ansatz 

.iCii,niYjj^{(^,4')j with Y^{9,(j)) a spherical harmonic, a set of linear 



equations for the coefficients ai^m is obtained by substituting (59) into (58) and expanding 
both sides in spherical harmonics. Solving the resulting linear equations always yields a 
valid P-representation. However, such P-representation constructions are non-unique in 
consequence of the identity [153| 1154] 



/ d'^x Y^{x) \x){x\ = for all integer I > 2j, 

Jin 



(60) 



which can be proved directly from (56) as a consequence of the addition law for angular 
momenta. This result shows explicitly that the set of coherent states \x) is over-complete 
(as is well-known). 

With the above as background, we now slow the pace of presentation. We consider 
the problem of finding a positive P-representation for a given operator p, that is to say, a 
representation for which P(x\p) > for all x. 

Positive P-representations have the useful property (for simulation purposes) of allowing 
us to interpret P(x\p) as a probability distribution over spin directions x. But from a 
mathematical point of view, distressingly little is known about positive P-representations, 
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in the sense that there is no known general method for constructing them, or even for 
determining whether they exist in a given case. 

We will now consider an operator that often appears in practical QMOR simulations: 
the thermal operator p^j^ defined by 

pf = exp(-/3t-s) , (61) 

where {si, S2, S3} are the usual spin-j operators satisfying [si, S2] = is^ (and cyclic permuta- 
tions), and i is a unit axis along which a spin is thermally polarized with inverse-temperature 



(3. By inserting a complete set of states into (58 ) and then substituting Wigner's expression 
(56), a well-known |153| 1162] closed-form expression for the Q-representation of p^j^ can be 



obtained: 

j 3 

Qmpf) = ^ ^ {x\j,m){j,m\pf\j,m'){i,m'\x) 

m=—j m'=—j 

= (coshi/3-a;-tsinhi/?)2^ (62) 

We now exhibit a closed-form analytic expression for a positive P-representation of the 
thermal operator as a distribution over coherent states, which is exact for all j and /?: 

Design Rule 3.2. A positive P -representation for the spin-j thermal operator pj^ is given 
in terms of the Q -representation by 

P{x\pf) = 1/Q{-x\pT+i) . 

To our knowledge this is the first such P-representation given, other than the j — > 00 limit 
of quantum optics in which P and Q are both simple Gaussians. 

If we regard the above result solely as a mathematical theorem to be proved by the most 
expedient means, we can do so by treating it as an ansatz. The resulting proof is short. 



Taking matrix elements of the defining relation (59) betwe en s tates {j,m\ and \j,m'), and 
without loss of generality setting t = (0, 0, 1), Design Rule 3.2 is equivalent to the definite 
integral 

2.?' + ^ f {j,rn\x){x\j,m') 



Substituting the Wigner representation (56) and the Q-representation (62) into (63), we 



can check by numerical integration that (63) is correct for randomly chosen values of j, m. 



m' , and /3. Thereby encouraged, we soon discover an integration strategy by which (63) 
yields to analytic evaluation in the general case. In brief, the 1?^ angle integration yields 
the requisite Smm' factor; the 0-angle integration can be transformed into an integral over 
rational functions in 2; = cos0 via identities like cos^-' ^9 = {1 -\- zy /2^; and the resulting 
integral is recognizably a representation of the Gauss hypergeometric series [U eqs. 15.1.8 



& 15.3.8] that evaluates to (63). 



Design Rule 3.2 is simultaneously frustrating, reassuring, and intriguing. It is frustrating 
because our QMOR analysis will answer the natural question "Where did the ansatz come 
from?" by "It is the solution to a Fokker-Planck equation that describes spin-systems 
in thermal equilibrium." But this answer provides no satisfying rationale for why the 
P-representation exists, or for its simple analytic form. Furthermore, QMOR analysis will 
provide no answer at all to the natural follow-on question, "Is there a reason why the thermal 
operator's positive P-representation is simply given in terms of its Q-representation?" 
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Design Rule 3.2 is reassuring because it tells us that QMOR simulations can simulate 
thermal equilibrium (at least in simple cases). From a physical point of view this reassures us 
that the dimensional reduction associated with QMOR mechanics preserves at least some 
crucial thermodynamic physics. From a practical point of view it provides a reassuring 
consistency check that the (rather lengthy) chain of theorems and stochastic analysis that 
leads to the P-representation is free of algebraic errors. 



And finally, Design Rule 3.2 is intriguing because it suggests that a partial answer to the 
open question "Which operators have positive P-representations?" might be "Only thoses 
operators that are proportional to a density matrix that is associated with a stationary 
measurement-and-control process that drives input states to coherent states." This conjec- 
ture is intriguing because its negation is logically equivalent to "Density matrices exist that 
cannot be the result of any stationary measurement-and-control process that drives input 
states to coherent states." Such density matrices would be remarkable entities from both 
a mathematical and physical point of view, and in particular their existence (or nonexis- 
tence) is directly relevant to the further development of practical design rules for QMOR 
simulations. The conjectured answer is therefore intriguing whether it is true or false. 

3.2.3. QMOR respects classical linear control theory QMOR simulations of macroscopic 
objects (like MRFM cantilevers) regard them as spin-j quantum objects having very large 
j. We will see that the resulting dynamics typically are linear. Engineers have their own 
idioms for describing linear dynamical systems, which are summarized in block diagrams 
like this one: 



/ext -^(y^ Go ^Or^ <r (64) 



H 



To show that QMOR simulations can accurately model systems like the above, we will trans- 
form the above diagram — using strictly classical methods — into an operationally equivalent 
form that more naturally maps onto QMOR algorithms and geometric structures. 

For convenience, these classical equivalences are summarized as design rules in Figs. [6]- 
[Tj Experienced researchers will recognize these equivalences as being elementary, but to 
the best of our knowledge, they have not previously been recognized in the literature of 
engineering or physics. 

As with the Feynman diagrams of physics, the block diagrams of engineering depict 
systems of equations. Our diagram conventions are standard, as briefly follows. The block 



diagram (64) corresponds to a set of linear relations between force noise f^{t), measurement 
noise (/"(t), input external force f^^^{t), and output measurement q™{t), which we take to 
be classical real-valued functions. In particular we specify that f^{t) and q^{t) are white 
noise processes having correlation functions C satisfying 

C{q''{t)q''{t')) = \ Sgn5{t - t') (65a) 
C{f-{t)f^{t')) = lSfn6{t-t') (65b) 
C{r{t)q^{t')) = (65c) 

so that Sqn and Sfn are (one-sided) white-noise spectral densities. A circle *p± is a node 
whose inputs are added and subtracted, a cross (S) is a node whose inputs are multiplied, a 
triangle indicates a positive real scalar gain 7, and a square box _K_ indicates convolution 
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CLASS 


BLOCK 


DEFINITION 


REMARKS 


dynamical: 




G{u) = G{-oj) 


' Dynamical kernels are 
^ time- reversal invariant. 


feedback: 




H[t) = for r < < 


' Feedback kernels are causal 
(including feedback from 
thermodynamic reservoirs). 


backaction: 




'Hioj) = i sgnu 
7 > is real 


Backaction kernels 
are Hilbert transforms 
followed by a gain 7 



Figure 6: The three kernel classes of linearized QMOR simulations. 



with a general real-valued stationary kernel K such that 

/oo 
dt' K{t-t')a{t') ^ a^VKWh (66) 

(is depicted as) 

Alternatively, convolution blocks can be specified in the Fourier domain. Our Fourier trans- 
form convention is that a{uj) is defined to be 



00 



a{uj) = / dr e~"^^a{t) (67) 



00 



and similarly for b{uj), K(uj), etc. Therefore a frequency-domain description of (66) is 



b{Lo) = K{uj)d{uj) a-^{KWb (68) 

(is depicted as) 

We build our QMOR block diagrams from three classes of linear classical kernels: dynam- 
ical kernels (-^[g^]-^), feedback kernels (-^[h]-^), and backaction kernels (-»-[?i]->-), whose 
defining mathematical properties are specified in Fig. [6j 

In brief, dynamical kernels by definition are time-reversal invariant, feedback kernels by 
definition are causal, and the backaction kernel is the Hilbert transform that is well-known 
(and much-used) by signal processing engineers. We remark that the Hilbert transform 
is formally non-causal, but in practical narrow-bandwidth applications (like radio trans- 
mitters, acoustic processors, MRFM cantilever controllers, etc.) its effects can be closely 
approximated by a causal derivative transform. 

For causal kernels, analytic continuation from the Fourier variable u; to the Laplace 
variable s = iu is well-defined. Partly because causal kernels are of central importance in 
control engineering, and partly by tradition, Laplace variables are more commonly adopted 
in the engineering literature than Fourier variables (although both are used). However, the 
Laplace analytic continuation of the non-causal Hilbert transform kernel Ti^uj) = isgn(uj) is 
not well-defined. For this reason our analysis will focus exclusively upon upon time-domain 
and frequency-domain (Fourier) kernel representations. 

Central to QMOR analysis and simulation is the purely mathematical fact (which also 
appears in Fig. [7] as Design Rule 3.4) that the following classical systems are operationally 
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equivalent: 



1 



H 



= r 




(69) 



1/2 



By "operationally equivalent" we mean that the dynamical relation between the applied 
forces f^^^{t) and the measured position g™(t) is identical for the two sets of equations, as 
are the stochastic properties of g™(t). The internal state of the system is of course very 
different for the two cases, but so long as we adhere to the strict operational principle 
"never observe the internal state of a system (even a classical system)" this difference is 
immaterial. 

Developing design rules for QMOR mechanics is considerably easier if we habituate 
ourselves to the Hilbert transform that appears in all three classical design rules of Fig. [7j 
We begin by remarking that from an abstract point of view, the Hilbert transform defines 
a complex structure on the space TZ of real- valued functions r, that is, a map TC: TZ ^ TZ 
that satisfies Tl^ = —I. This corresponds to two diagrammatic identities 



n-n 



and 



n 



rv 



n 



(70) 



that make it obvious (since the right-hand sides are causal) that the non-causal aspects of 
-^[^]-^ vanish when it is applied an even number of times. Because the even-numbered 
moments wholly determine the statistical properties of (zero-mean) Gaussian noise, we 
begin to see how it is that noncausal Hilbert transforms can appear in the equivalences of 
our classical design rules without inducing observable causality violations. In the broader 
context of quantum simulations, causality is assured by the Theorema Dilectum, and we 
will establish in a later section that the Theorema Dilectum is directly responsible for the 
appearance of Hilbert transforms in the classical limit. 



We saw already (in Section 1.5.3) that complex structures are a defining geometric 



characteristic of Kahler manifolds, so the appearance of Hilbert transforms in real- valued 
classical dynamics is a mathematical hint that noisy classical dynamical trajectories support 
a complex structure that projects naturally onto gabion-Kahler manifolds. This complex 
structure is manifest not in the (real-valued) state variables of classical systems, but in the 
causal properties of the response of classical systems to noise. 

3.2.4- Remarks upon the spooky mysteries of classical physics For teaching purposes, it is 
helpful (and amusing) to pretend that we live in a world in which linear control theory is 
taught according to a nonstandard ontology in which the Hilbert transform has a central 
role. This ontology was conceived as a philosophical provocation (a herausforderung |175j ). 
but it has subsequently proved to be useful in teaching and a fertile source of technical 
inspiration. We will call it the classical Hilbert ontology, or sometimes just the Hilbert on- 
tology. Our motivation for emphasizing the mysterious properties of classical measurement 
theory is similar to the motivation of Wheeler and Feynman |196| 1197] in proposing their 
non-causal classical electrodynamic ontology proposed in the 1940s. 
The tenets of Hilbert ontology are taken to be: 
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Three design rules that reflect the "spooky mysteries of classical physics" 

Design Rule 3.3. The net action of two independent white noise processes that act upon 
opposite sides of a linear dynamical kernel -»-[G^i]->- is operationally equivalent to a single white 
noise process that acts with (noncausal) Hilbert hackaction- and- gain --\h\-^>—' as follows: 




Design Rule 3.4. A linear dynamical system that is subject to independent force noise and 
measurement noise is operationally equivalent to the same linear dynamical system acted upon 
only by measurement noise having Hilbert backaction as follows: 




Design Rule 3.5. In narrowband dynamical systems, defined as those having a Hilbert 
gain 7 satisfing 7 ^ /c, where k is the scale of the system spring constant, it is feasible in 
practice to operationally approximate zero-noise internal system dynamics as follows: 

to—?" 



r 



H 



^O- Go 



H 



1/2 

n 



7 = 



1/2 



(for in-band frequencies) 



Figure 7: Three design rules that reflect "spooky mysteries of classical physics." 

The above design rules follow from elementary classical considerations as (briefly) follows. 
Design Rule 3.3 follows from the time-reversal invariance of the dynamical kernel Go (see 
Fig. [g]) plus the statistical properties of Gaussian noise. Design Rule 3.4 is then obtained 
by adjoining a causal feedback kernel H to the diagram of Rule 3.3. Design Rule 3.5 then 
follows as a practical application of Rule 3.4, in which a causal kernel H approximates a 
non-causal Hilbert backaction kernel Ti. within a device's finite operational bandwidth. 
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• By Occam's Razor, ontologies that invoke one noise source are preferred over ontolo- 
gies that invoke two. Therefore the Hilbert ontology regards backaction as being a 
physically real phenomenon that is universally present in all noisy dynamical systems 
(its mathematical expression being Design Rule 3.3 of Fig. [?]). 

• It follows (by Design Rule 3.4) that measurement noise always back-acts upon system 
dynamics in such a way as to "drag" the state of the system into agreement with the 
measurement. This state-dragging Hilbert backaction has a central ontological role: it 
is the fundamental mechanism of nature that causes measurement processes to agree 
with reality. 

• As a measurement process approaches the zero-noise limit, the increasingly strong 
state-dragging Hilbert backaction from that measurement process dynamically "col- 
lapses" the variable being measured, so as to force exact agreement with the measure- 
ment. 

• As a corollary, zero-noise measurements are unphysical, because they imply infinitely 
strong backaction. That is the explanation in Hilbert ontology of why all real-world 
measurement processes are noisy. 

• In narrow-band systems, it is possible to cancel the Hilbert backaction noise via causal 
feedback control. This means that zero-temperature narrow-band systems can be 
simulated, or even realized in practice, provided that all noise processes are accessible 
for purposes of feedback (as mathematically expressed by Design Rule 3.5 of Fig. [7|. 

• Although the mathematical kernel associated with Hilbert feedback is non-causal, 
this noncausality cannot be exploited for purposes of communication, for the physi- 
cal reason that the backaction kernel transmits only noise. The mathematical proof 
is simply that Design Rules 3.3-5 can be expressed as equivalent systems having all 
kernels causal (as shown in Fig. [7]). But in Hilbert ontology this causality proof is 
viewed as being a purely mathematical artifice, because it postulates "spooky" un- 
correlated forces that "obviously" have no physical reality (according to our Occam's 
Razor tenet) since no mechanism is given for them. 

For purposes of this article, we designate the above Hilbert ontology to be the "true" classi- 
cal reality of the world, and we seek to provide a microscopic justification of it from orthodox 
quantum mechanics (recognizing of course that the Hilbert ontology has been constructed 
specifically to ensure that our overall herausforderung yields practical and interesting re- 
sults). 

From a teaching point of view, the Hilbert ontology helps students appreciate that the 
mysteries and ambiguities that traditionally are taught as belonging exclusively to quantum 
mechanics — like "wave functions collapsing to agree with measurement" and "noncausal 
correlations" — are manifest too in (at least one) wholly classical ontology. 

It is traditional in both the popular and the scientific literature to call these ontological 
mysteries and ambiguities "spooky." The popularity of "spooky" in the scientific literature 
can be traced to an influential article by Mermin |138j . who adopted it from an idiomatic 
phrase "spukhafte Fernwirkungen" that Einstein uses in the Einstein-Bohr correspondence; 
this phrase is generally translated as "spooky action at a distance" [321 11U51 1141| . Given 
the central importance of the spookiness of quantum mechanics, it would be astonishing if 
this spookiness was wholly invisible at the classical level. 

Our Hilbert ontology therefore serves both to guide our calculations and to help us 
celebrate the "spooky mysteries of classical physics." 
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3.2.5. Experimental protocols for measuring the Hilbert parameters A crucial test of an 
ontology is whether it motivates us to proceed to practical calculations that yield useful 
and/or surprising results; we now do so. We consider a laboratory course in which students 
are guided by the Hilbert ontology to explore the fundamental and practical limits of low- 
noise sensing and amplification. For definiteness, we assume that the students work with 
nanomechanical oscillators as force detectors (as in MRFM technology), but a similar course 
could feasibly be organized around radio-frequency (RF) sensing and amplification, optical 
sensing and amplification, or acoustic sensing and amplification. 

We consider the experimental problem of measuring the two Hilbert parameters of a 
nanomechanical oscillator, namely the measurement noise Sq and the Hilbert gain 7. We ne- 
glect the intrinsic damping of the oscillator, such that the dynamical kernel is Go(i-i;) = 
l/[m{LL!Q — tj^)], where m is the mass of the oscillator and luq is the resonant frequency. 

The measurement is readily accomplished by the following protocol. Derivative feedback 
control is applied having a kernel H(ij) = iTuj, and the value of the controller gain F is 
adjusted until the spectral density Sqm of the measured cantilever displacement is 
observed to be flat in the neighborhood of the cantilever frequency loq. The control gain 
adjustment is straightforward: if the spectrum of q"^{t) has a peak, then F is too small; if 
the spectrum has a hole, then F is too large. 

The preceding protocol is perfectly feasible in practice {i.e., it is not a Gedankenex- 
periment). The protocol fails only when the required dashpot gain is impractically large, 
such that Tljq > k, where k = moj^ is the spring constant of the cantilever. Such failures 
typically indicate that a measurement process has a noise level that is too large to be of 
practical interest. 

In practical cases the required control gain typically satisfies Fq <C F <C A;, where Fq = 
k/{ujQQ) is the intrinsic damping of a cantilever having quality Q. In such cases the intrinsic 
damping Fq has negligible effect on the time-scale of the controlled response, and so we can 
regard the cantilever's dynamical kernel as having the time-reversal-invariant dynamical 
form Go described in Fig. |6] The design rules of our Hilbert classical ontology therefore 
can be applied without modification. Specifically, Design Rule 3.5 determines the Hilbert 
backaction to be 7 = Fwq, and determines the measurement noise to be 5gn = S'gm(a;o). 

Of course. Design Rule 3.5 mathematically assures us that the adherents of traditional 
classical ontology can — without operational contradiction — ascribe these same measure- 
ments to a fictitious force noise /"(t) having spectral density S'jn = 7^5qn, so that the 
vexing question of whether this force noise is "fictitious" or "real" is operationally immate- 
rial. 

3.2.6. QMOR simulations respect the fundamental quantum limits We now apply the re- 
sults of the preceding section to establish criteria that QMOR simulations must satisfy 
in order to respect the known fundamental quantum limits to amplifier noise and force 
measurement. 

We suppose that a force signal f{t) = /o cos(a;o^ + <Po) is applied to the oscillator. The 
carrier frequency ujq is tuned to match the oscillator resonance frequency, and the unknown 
magnitude /o and unknown phase (f)Q of the signal are to be determined from measurement. 
This is a common task in practice. 

When the oscillator is configured according to the measurement protocol of the preceding 
section, then according to the right-hand block diagram of Design Rule 3.5 of Fig. [7] the 
mean power po absorbed by the oscillator's (wholly classical) feedback controller during 
the measurement process is po = /oV(27). The absorbed power inferred from the (wholly 
classical) measurement record (?™(t) has an equivalent (one-sided) noise PSD S'p"*''"* whose 
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expression in terms of the Hilbert parameters is readily shown to be 

= 4poLOo7S,n (71) 

Now we connect this result with a known fundamental quantum limit on noise in power 
amplifiers. We adopt the definition of Caves [45^: 

An amplifier is any device that takes an input signal, carried by a collection of bosonic 
modes, and processes the output to produce an output signal, also carried by a (possibly 
different) collection of bosonic modes. A linear amplifier is an amplifier whose output 
signal is linearly related to its input signal. 

Following a line of reasoning put forth in the 1960s by Heffner |106] and by Haus |104] . 
which Caves article develops in detail |l5] , we suppose that the input power is supplied by 
a bosonic mode- type device (like a resonant circuit, or an RF wave-guide, or a single-mode 
optical fiber) whose power level has been independently measured with a shot-noise-limited 
photon counting device. According to orthodox quantum mechanics such counting processes 



have Poisson statistics (we will establish in Section 3.4 that QMOR simulations respect this 
rule), and therefore the input power has a quantum-limited noise PSD given by 

S^/"* = 2^,0^0 (72) 

The measured power and phase suffice to create an arbitrary-gain replica of the input signal, 
and the noise- figure (NF) of this effectively infinite-gain power amplifier is simply given from 



(7lj) and ^ by 

NF = 57*P^V4"P"* = 2jSqn/h = 2(5/n5gn)V2/^ (73) 

In Caves' nomenclature, we are regarding our continuously measured oscillator as an equiva- 
lent "phase-insensitive linear amplifier" having infinite gain. The analyses of Heffner, Haus, 
and Caves |45t I1U4| 1106] establish that what Caves calls the "fundamental theorem for 
phase-insensitive power amplifiers" is simply NF > 2 (in the infinite-gain limit), which in 
decibels is lOlog^Q 2 ~ 3 dB. 

We remark that the 3 dB quantum limit to power amplifier noise has been experimentally 
observed |108j and theoretically analyzed |176] since the early days of maser amplifiers in the 
1950s; our review focuses upon the work of Heffner, Haus, and Caves solely because their 
analyses are notably rigorous, general, clearly stated, and (importantly) their predicted 
quantum limits are mutually consistent and consonant with subsequent experiments. 



Our result (73) then establishes that the Heffner-Haus- Caves noise-figure limit finds its 
expression in QMOR analysis in the following three equivalent ways: 

NF > 2 ^ SfnSgn >h^ ^ -fSgn > H (74) 

THE HeFFNER-HAUS-CaVES LIMIT THE BrAGINSKY-KhALILI LIMIT THE HiLBERT LIMIT 

The middle expression we recognize as the continuous-measurement version of the standard 
quantum limit to force and position measurement, in precisely the quantitative form derived 
by Braginsky and Khalili , which in turn derives from earlier seminal work by Braginsky, 
Vorontsov, and Thorne [23]. Although the Braginsky-Khalili limit was derived by very 
different methods from the Heffner-Haus-Caves' limit, we see that the two quantum limits 
are equivalent. To the best of our knowledge, this equivalence has not previously been 
stated in the above quantitative form. What we have chosen to call "the Hilbert limit" on 



the right-hand side of (74) has (to our knowledge) not previously been recognized anywhere 
in the literature, and yet from the viewpoint of Hilbert ontology it is the most fundamental 
of the three. 

We express the above three-fold equivalence as a fundamental QMOR design rule: 
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Design Rule 3.6. QMOR simulations respect the quantum measurement limit in all its 
equivalent forms: the noise-figure (NF) limit in power amplifiers (NF > 2), the standard 
quantum limit to the measurement of canonically conjugate variables (SfnSqn > f^), and 
the Hilbert limit that measurement noise and state- dragging Hilbert backaction cannot both 
be small ('^Sqn > h). 



We present a quantum derivation of these limits in Section 3.4, in the context of a more 
general analysis that encompasses nonlinear quantum systems. 

3.2.7. Teaching the ontological ambiguity of classical measurement What should we teach 
students about the internal state of the system during the preceding (strictly classical) 
protocols for measuring the Hilbert parameters Sqn and 7? 

From a strictly logical point of view, this question need not be answered, since our 
measurement protocols and design rules are careful to make no reference to the internal 
state (even at the classical level). But in practice, some answer must be given, for the 
pragmatic reason that students require some coherent story about what the systems they 
are measuring are "really" doing. This coherence is especially necessary to those science 
and engineering students (a substantial portion, in our view) whose careers will require that 
they extend their professional expertise from the classical to the quantum domain. 

It is therefore desirable that common-sense questions from students receive common- 
sense answers from teachers, and it is equally desirable that these answers prepare for a 
classical-to-quantum educational transition that is as nearly seamless as is practicable. 

Adherents of traditional classical ontology can argue cogently as follows: "It may be op- 
erationally correct to ascribe experimental results wholly to measurement noise and Hilbert 
backaction, but it is physically wrong, because we have strong physical reasons to believe 
that the excitation of the oscillator is being driven by force noise from a thermal reservoir 
whose internal dynamics we do not observe." 

Adherents of Hilbert classical ontology can offer an similarly cogent counter- argument, 
which however requires the assertion of a definite mathematical result (in italics): "We too 
believe that the observed excitation of the oscillator is in fact being driven an unobserved 
thermal reservoir, but the action of a unobserved thermal reservoir can be ascribed to a 
covert process of measurement, Hilbert backaction, and control." 

The practical consequence of the above reasoning is that this article's embrace of Hilbert 
ontology has practical utility only if we can develop well-posed QMOR algorithms by which 
the action of thermal reservoirs upon a dynamical system is simulated by unobserved / covert 
processes of quantum measurement, Hilbert backaction, and classical control. We develop 



the necessary algorithms in Section 3.4 using as our main mathematical tool the Theorema 



Dilectum that was already given as Design Rule 3.1 A key mathematical result will be the 
positive P-representation that was already given as Design Rule |3.2[ this helps us foresee 
the analysis path by which the QMOR analysis program will succeed. And of course Design 
Rules 3.3-5 will emerge naturally in the course of our analysis. 

In summary, for students and teachers alike, a sufficient justification for embracing the 
Hilbert ontology is that it leads to QMOR simulation algorithms that are computationally 
efficient, operationally orthodox, and mathematically novel. A provocative side-effect is 
that we learn to perceive the "spooky mysteries of quantum physics" as being manifest in 
all noisy dynamical systems, even classical ones. 

3.3. Physical aspects of QMOR 

We now turn our attention to the physical aspects of measurement, our goal being to 
establish connections between the concrete description of measurement in terms of hardware 



3.3 Physical aspects of QMOR 



55 



and experimental protocols on the one hand, and the measurement operator formalism of 
our QMOR algorithms on the other hand. 

3.3.1. Measurement modeled as scattering We adopt as the fundamental building block of 
our simulations a single particle of spin j, described by a wave function I'i/^), and numerically 
encoded as a complex vector with dim \ = 2j + 1. We will simulate all noise and all 
measurement processes by scattering photons off the spin, one photon at a time, and we 
associate each scattering event with a single time step in Fig. [T] 

We envision photon scattering as the sole mechanism by which noise is injected into 
our simulations, and interferometry as the sole means of measurement. We describe photon 
scattering as a unitary transformation acting on the spin state 

(before scattering) \ipn) — > exp(i20s°P) (after scattering). (75) 

A purely conventional factor of 2 is inserted in the above to simplify our calibration rules. 

We will call s°p a measurement generator. In general s°p can be any Hermitian matrix, 
but in our simulations it will suffice to confine our attention to s°p E {si, 52,83}, where 
{si, S2, S3} are rotation matrices satisfying the commutation relation [sj, s^] = iejkisi. These 
matrices generate the rotation group, and our discussion will assume a basic knowledge of 
their algebraic properties, which are discussed at length in many textbooks {e.g., |93| I139| 

nmnrm\ ). 

We adopt the near-universal convention of working exclusively with spin operator matrix 
representations that are irreducible and sparse, having dimension 2j + 1 for j the spin 
quantum number, with S3 a real diagonal matrix, si real and bidiagonal, and S2 imaginary 
and bidiagonal. The scattering strength is set by the real number 6, which in our simulations 
will always satisfy <C 1. 

3.3.2. Physical and mathematical descriptions of interferometry Figure [s] provides three 
idealized diagrams for thinking about the physical, geometric, and algebraic aspects of 
photon interferometry, similar to the idealized diagrams in the Feynman Lectures |69| chs. 5- 
6] that depict the Stern-Gerlach effect. 

Figure [sj^a) depicts physical fiber-optic interferometers that confine photons within low- 
loss optical fibers |26| 1171] . The region of overlap between the fibers allows photons to 
transfer from one fiber to another with an amplitude that is subject to engineering control. 
This overlap region plays the role of the semi-silvered mirror in a traditional Michelson 
interferometer. 

Figure [sj^b) depicts optical couplers geometrically, as general-purpose devices for link- 
ing incoming and outgoing optical amplitudes. Couplers thus can be braided into optical 
networks of essentially arbitrary topology. 

Figure [sj^c) depicts optical couplers algebraically, as linear maps between incoming com- 
plex optical amplitudes a™ = {aj.", at?) ^S' ^brl outgoing optical amplitudes a°^* = 
{al^\ a'^f-, agf , a^f } such that a°"* = S • a^"^, with S the optical scattering matrix. With 
regard to Fig. [8]^c), the index "tl" is the "top left" port, "br" is the "bottom right" port, 
etc. . Our normalization convention is such that the probability of detection of an outgoing 
photon is |a°"*p. Optical losses in real- world couplers are small, such that ja"*^*!^ = |a™P, 
i.e. the optical scattering matrix S is unitary. 

3.3.3. Survey of interferometric measurement methods The stochastic measurement and 
noise processes in our simulations can be conceived as interferometric measurements per- 
formed on each scattered photon, and this physical picture will prove very useful in designing 
MOR techniques. Such measurements require that each incoming photon be interferometri- 
cally split before it scatters from the spin, to allow subsequent interferometric recombination 
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Figure 8: Three aspects of photon interferometry. 

Photon interferometry from the (a) physical, (b) geometric, and (c) alge- 



braic points of view. See Section 3.3.2 for discussion 



and measurement. It is convenient to conceive of this initial splitting as performed by a 2 x 2 
single-mode optical fiber coupler, as illustrated in Figs. |9]^a-c). The devices of this figure 
may be regarded as physical embodiments of the simulation algorithm of Fig. [T] 

This physical picture embodies the idealizing assumptions that optical couplers are 
exactly unitary, that photon emission into the fiber takes place at equally spaced intervals 
6t, and that photon detectors register a single classical "click" upon detection of each 
photon, which occurs with detection probability |a°"*p. The unitarity of photon scattering 
and interferometric propagation then ensures that each incoming photon results in precisely 
one detection click. 

With respect to the algorithm of Fig. [ij the stochastic selection of operator M^^j versus 
M(^j is physically identified with these clicks, such that the sole data set resulting from 
a simulation is a set of classical binary data streams, with each stream comprising the 
recorded clicks for a measurement operator pair. Such binary streams closely resemble real- 
world MRFM experimental records, in which signals corresponding to photoelectrons from 
an optically-monitored cantilever are low-pass filtered and recorded. 

We remark that the interferometers of Fig. [9]^a-c) can be regarded with equal validity 
either as idealized abstractions or as schematic descriptions of real- wo rid experiments. For 
example, in our simulations we will regard MRFM cantilevers as spins of large quantum 
number j, in which case the interferometer geometry of Fig. [oj^c) is identical (in the topol- 
ogy of its light path) to the real-world fiber-optic interferometers used in typical MRFM 
experiments [26l fTTnl ITTT] . 

Trapped-ion experiments are examples of continuously observed small-j quantum sys- 
tems, since the quantum mechanics of a two-state atom is identical to the quantum mechan- 
ics of a spin-^ particle. Such experiments presently are conducted with photons detected 
directly as in Fig. |18| I130| [T58]. We remark that it would be theoretically interesting, 
and experimentally feasible, to conduct trapped ion experiments with weakly interacting 
photons detected interferometrically, as in Fig. |9|^c). 

In trapped ion experiments, the observed transitions in the photon detection rates are 
observed to have random telegraph statistics, which are typically attributed to "quantum 
jumps" |18| I130j or to "instantaneous transitions between energy levels" jl58l . As was rec- 
ognized by Kraus more than twenty years ago |123| p. 98], "such reasoning is unfounded," 
and we will see explicitly that observation of telegraph statistics does not imply discontin- 
uous evolution of \ip). Thus, readers having a classical MOR background need not regard 
discussions of quantum jumps in the physics literature as being literally true. 
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Figure 9: A physical illustration of the Theorema Dilectum. 

This physical embodiment of the formal QMOR simulation algorithm of Fig. [T] using 2x2 
fiber-optic couplers to scatter photons off a spin state. In (a) the photons are detected 
with equal probability, such that no information about the state is obtained. In (b) the 
phase is detected via homodyne (self-interfering) interferometry, such that a small amount 
of information about about the state is obtained. To an outside observer, the fate of the 
downstream photons is immaterial, hence (a) and (b) must embody physically indistin- 
guishable noise processes, despite their differing quantum description. This is the physical 
content of the Theorema Dilectum. From a mathematical point of view, the free choice 
of an arbitrary downstream unitary transform upon the photon paths is manifested in the 
[/-transform invariance of (47). (On-line version only: photons in red are causally down- 



stream of the scattering, and thus can be measured in arbitrary fashion without altering 
the physically observable properties of the noise being simulated.) 



3.3 Physical aspects of QMOR 



58 



3.3.4- Physical calibration of scattering amplitudes The calibration of our simulations will 
rely upon a physical principle that is well-established, but somewhat counterintuitive. The 



principle is this: measurements of the scattering phase (75) suffice to provide detailed infor- 
mation about the Hamiltonian that is responsible for the scattering. Gottfried's discussion 
of Fredholm theory in atomic scattering [931 Section 49.2] provides a good introduction to 
this topic. 

As a calibration example, we consider here the measurement process associated with 
the simple interferometer of Fig. |9]^a). Using the S-matrix of Fig. |8]^c) to propagate the 
input photon in Fig. 9|a) through the apparatus, the measurement operators {Adj^^^ , Alj^^^} 
associated with detection on the (+) and (-) channels are 
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(76) 



The above equation follows the convention that the optical amplitudes of the top (bottom) 
fiber path in Fig. |9](a) are listed as the top (bottom) element of the above column arrays, 
such that successive interferometric couplings are described by successive 2x2 unitary ma- 
trix operations. This convention is common in optical engineering because it unites the 
geometric and algebraic descriptions of Fig. [sj^b-c). In the context of our quantum sim- 
ulations, each element of the above arrays is formally an operator on the wave function 
IV'), but for c-number (complex number) array elements like "0", "1" and "z" an implicit 
identity operator is omitted for compactness of notation. These identity operators physi- 
cally correspond to events that are not dynamically coupled to the spin state, like photon 
emission, propagation through interferometers, and subsequent detection. 

The resulting overall measurement operators {m\^^ , M^^^^ describe the state-dependence 
of the scattered phase of the detected photon, and it is evident that they satisfy the mea- 



surement operator completeness relation ( 45 ) . 



3.3.5. Noise- induced Stark shifts and renormalization We now consider an optical scatter- 
ing effect known as the AC Stark shift, as induced by the scattering process of Fig. [9]^a). 
Applying the simulation algorithm, we find that during a time At ^> 6t the final state of 
the simulation accumulates a state-dependent phase shift such that 



|V(At))=exp(i2n„es°P)|V'(0)), 



(77) 



where n(_,.) is the number of photons detected on the (+)-channel. It is easy to show that 
n{+) has mean = At/{26t) and standard deviation = y^/x^!,!^. 

We see that the average effect of the photon scattering is equivalent to a dynamical 
Hamiltonian H°^, such that 

H°P = -es°Ph/6t, (78) 

which we identify as the effective Hamiltonian of a Stark shift. The Stark shift fluctuates 
due to statistical fluctuations in the number of photons detected on the (+)-channel, such 
that in the continuum limit the photon detection rate r(t) is a stochastic process having 
mean = l/{26t) and white-noise spectral density Sr = l/(45t) (in a two-sided spectral 
density convention). This implies that the Stark shift fluctuations have an operator- valued 
power spectral density (PSD) 



Shop = e^s^Pfh^/dt = {H°Pf5t. 



(79) 



This result calibrates the externally-observable Stark shift parameters {H°p, Sh°-p} in terms 
of the internal simulation parameters {9s°P,5t} and vice versa. 



3.3 Physical aspects of QMOR 



59 



The preceding two equations (78-79) reflect the weh-known phenomenon in physics 
that interaction with a measurement process or thermal reservoir renormalizes the physical 
properties of a system. But as presented here, these same equations exhibit a known 



pathology in the limit 6t — > 0: if we take 6 cc V St such that the Stark noise (79 ) has a finite 



limit, then the magnitude of the Stark shift Hamiltonian (78) diverges. The origin of this 
(unphysical) infinite energy shift is our (equally unphysical) assumption that measurement 
"clicks" occur at infinitely short intervals, such that the PSD of the measurement noise is 
white; such white noise processes inherently are associated with infinite energy densities. 

In the present article we will simply repair this white-noise pathology by simply adding 
a counter-term to the measurement process, such that the Stark shift is zero in all our 
measurement processes. In the language of renormalization theory, we redefine all of our 
measurement operators so that they refer to the "dressed" states of the system. See {e.g.) 
|180] and |128j for further discussion of this delicate point, whose detailed analysis is beyond 
the scope of the present article. 



3.2.1 



we discussed the Theorema 



3.3.6. Causality and the Theorema Dilectum In Section 
Dilectum from a mathematical point of view that defined it terms of the (unitary) U- 
transformation of (47). In Fig. |9]^a-c) we associate the [/-transformation with the physical 



choice of what to do "downstream" with photons that have scattered off the spin. If 
we regard this downstream choice as being made made by Alice, and we assume that 
Bob is independently monitoring the same spin as Alice, then the physical content of the 
Theorema Dilectum is that Alice's downstream measurement choices can have no observable 
consequences for Bob, and in particular, Alice cannot establish a communication channel 
with Bob via her measurement choices. 



We therefore physically identify U in (47) with the downstream optical couplers in 



Fig. |9|^a-c), such that the algebraic freedom to specify an arbitrary unitary transform U is 
identified with the physical freedom to "tune" the interferometer by adjusting its coupling 
ratio and fiber lengths (which adjust the phase of the output amplitudes relative to the 
input amplitudes). Because experimentalists are familiar with this process and with the 
phrase "interferometer tuning" to describe it, we adopt the word "tuning" to describe the 
process of adapting U to optimize {M(+), M(_)} for our simulation purposes. 

3.3.7. The Theorema Dilectum in the literature The invariance associated with U has 
received most of its attention from the physics community only recently; it is not mentioned 
in most older quantum physics texts. Both Preskill's class notes |159j and Nielsen and 
Chuang's text |146] give physics-oriented proofs of the Theorema Dilectum. The crux of 
all such proofs is that the choice of U does not alter the density matrix associated with 
the ensemble average of all possible trajectories, such that the choice has no observable 
consequences. 

Carmichael |40| p. 122] seems to have been among the first to use the now-widely-used 
term unraveling in describing quantum trajectory simulations, and he explicitly recognized 
that this unraveling is not unique: 

We will refer to the quantum trajectories as an unravelling of the source dynam- 
ics since it is an unraveling of the many tangled paths that the master equation 
evolves forward in time as a single package. It is clear [. . . ] that the unraveling 
is not unique. 

(emphasis in the original). Diverse points of view regarding the ambiguity of trajectory 
unraveling can be found in the physics literature. Rigo and Gisin |166| argue that it is 
central to our understanding of the emergence of the classical world, and they make their 
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case by presenting four different unravelings of a single physical process; we will adopt 
a similar multiple-unraveling approach in analyzing the IBM single-spin experiment (see 



Section 4.2 ). Percival's text adopts the equally valid but sharply contrasting view that |152| 
p. 46] : "the ambiguity [. . . ] is a nuisance, so it is helpful to adopt a convention which reduces 
this choice." Breuer and Petruccione [25] simply state a result equivalent to the Theorema 
Dilectum, without further comment or attribution. 

Preskill's course notes and Nielsen and Chuang's text both take a middle point of view. 
They briefly describe the Theorema Dilectum as a "surprising ambiguity" |159| sec. 3.3] that 
is "surprisingly useful" \14:6\ sec. 8.2.4]. The word "surprising" invites readers to think for 
themselves about unravelling, and the word "ambiguity" suggests (correctly) that the im- 
plications of this invariance are not fully understood. Nielsen and Chuang's text notes that 
up to the present time, the main practical application of the Theorema Dilectum has come 
in the theory of quantum error correction, where the freedom to choose unravelings that fa- 
cilitate the design of error correction algorithms has been "crucial to a good understanding 
of quantum error correction" |146l sec. 8.2.4]. 

In the context of quantum computing theory, Buhrman et al. |27] have exploited the 
informatic invariance of the Theorema Dilectum to show that certain logical gates that 
are essential to universal quantum computing, when made noisy, can be indistinguishably 
replaced by randomly-selected gates from a restricted set of gates that can be simulated 
classically. Their proofs build upon earlier work by many authors |3l [231 llOlj and in 
particular upon the Gottesman-Knill theorem f92]. To our knowledge, this is the first 
formal quantum informatic proof that additional noise makes systems easier to simulate. 

3.4- Designs for spinometers 

We now present some basic designs for measurement operators constructed from the fun- 
damental set of spin operators {si, S2, S3} for particles of arbitrary spin j. 

We call this family of measurement operators spinometers, and we will characterize 



their properties in such sufficient detail that in Section 4.2 we will be able to simulate the 
single-spin MRFM experiment both numerically and by closed-form analysis. Most of these 
results have not appeared previously in the literature. 

3. 4-1- Spinometer tuning options: ergodic, synoptic, and hatrachian Starting with the 



fundamental spinometer pair {M,J)\M(5}that was defined in d76| and physically illustrated 
in Fig. [oj^a), ergodic spinometers are constructed by applying the following tuning 
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(80) 



Note that immediately following photon detection, and independent of the channel on which 
the photon is detected, a compensating Hamiltonian 9s°^ is applied to cancel the mean 



Stark shift that was noted in Section 3.3.4 The fluctuating portion of the Stark shift is 
not thereby cancelled, and it follows that ergodic spinometers are well-suited for simulating 
physically realistic noise, e.g., random magnetic fields that decohere spin states. 

We construct hatrachian spinometers from ergodic spinometers by adding a downstream 
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coupler, as shown in Fig. ^b). Our phase and tuning conventions are: 



M. 
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(81) 



The upper- right output port is tuned to be as dark as possible, such that detection clicks 
occur only sporadically; in experimental interferometery this is called dark port tuning. 
Each click that a dark port records is accompanied by a discrete jump in the wave function, 
hence the name "batrachian" for these measurement operators. We will see that batrachian 
tuning is well-suited to the analysis of data statistics. 

We construct synoptic spinometers similarly, but with a different phase tuning, as shown 
in Fig. ^c). Algebraically our tuning convention is 
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(82) 



We will see that synoptic spinometers do provide information about the quantum state — 
hence the name "synoptic" — and also that they compress quantum trajectories. 

We can now discern the general strategy of QMOR analysis: we model physical noise in 
terms of ergodic operators, we predict data statistics by the analysis of batrachian operators, 
and we compress simulated trajectories by applying synoptic operators. 

3.4-2. Spinometers as agents of trajectory compression The following derivations assume 
a knowledge of basic quantum mechanics at the level of Chapters 2 and 8 of Nielsen and 
Chuang ^146j (an alternative text is Griffiths [95]), knowledge of coherent spin states at the 
level of Perelomov |153| eqs. 4.3.21-45] (alternatively see Klauder |118] or del Castillo [H5]). 
and knowledge of stochastic differential equations at the level of Gardiner |841 sec. 4.3] 
(alternatively see Rogers [168!). 

Generally speaking, the design rules of this section were first heuristically suggested by 



the Hilbert ontology of Section 3.2 then confirmed by numerical experiments, and finally 
proved by analysis. Some of the lengthier proofs would have been difficult to discover 
otherwise; this shows the utility of the Hilbert ontology backed-up by numerical exploration. 

Nowhere in the derivations of this section will we make any assumption about the 
dimensionality of the Hilbert space in which the trajectory {IV'n)} resides. Therefore we are 
free to regard |^n) as describing a spin-j particle that is embedded in a larger multi-spin 
Hilbert space. Thus all the theorems and calibrations that we will derive will be applicable 



both to the single-spin MRFM Hilbert space of Section 4.2 and to the large-dimension 



"spin-dust" spaces that we will discuss in Section [43 

3.4.3. Spinometers that einselect eigenstates We define a uniaxial spinometer to be a mea- 
surement process associated with a single pair of measurement operators having generator 
We can regard as an arbitrary Hermitian matrix, since in a uniaxial measure- 
ment there are no other operators for it to commute with. We consider ergodic, synop- 



tic, and batrachian tunings as defined in (80-81). Without loss of generality we assume 
tr(s°P)2 = dim]?/;; 
of the coupling 9. 



I.e.. 



the mean-square eigenvalues of s°p are unity, which sets the scale 
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For a general state \ipn) and general Hermitian operator we define the operator 
variance A„(s°p) to be 

A„(S°P) = {^n\ - (^n|s°P|Vn))' IV'n) (83) 

remarking that in a finite-j Hilbert space 

A fs°Pl I ^ ^ '^"^ eigenstate of s°p, 
" \ > otherwise. 

Physically speaking, the smaller the variance A„(s°p), the smaller the quantum fluctuations 
in the expectation value {'ijjn\s"'^\ipn) ■ We will now calculate the rate at which measurement 
operators act to minimize this variance. 

Considering an ensemble of simulation trajectories, we define the ensemble-averaged 
variance at the n-th simulation step to be E[An{s°^)]. The algorithm of Fig. [T] evolves this 
mean variance according to 



MA„+i(.°P)] = ^ ^ E 
k=i ie{(+),(-)} 



1 ('il^n\(Mhh°PM'y\'iljr, 
(^„|(Mf)t(s°P)2M^'|^„). ^"^"'^ 



(V^„|(M'=)tMf|V^„) 



(85) 



For compactness we write the increment of the variance as 5A„(s°p) = A„+i(s°p)] — A„(s°p). 
Then for ergodic, synoptic, and bratrachian tunings the mean increment is 

ergodic tuning, 

^[5A„(s°P)] = { -402£;[A2(s°P)] synoptic tuning, (86) 
—B"^ E[J^n{s'^^)\ batrachian tuning. 



These results are obtained by substituting in (85) the spinometer tunings of (80-81), then 



expanding in 9 to second order. Here J- is the non-negative function 

^ , op, _ ((V'n|(g°P)^|V'n) - (V'n|g°P|V'n)(V'n|(^°P)^|^n))' 
^" (^n|(.°P)2|V'„) • 



Each term in the sequence {i?[Ai], £'[A2], . . .} is nonnegative by (83), and yet for synoptic 
and batrachian tuning the successive terms in the sequence are non-increasing (because in 
(86) the quantities A^(s°p) and .7>i(s°P) are nonnegative and there is an overall minus sign 



acting on them); the sequence therefore has a limit. For synoptic tuning the limiting states 
are evidently such that A„(s°p) — > 0, while for batrachian tunings .7-"„(s°p) — > 0, which in 
both cases implies that the limiting states are eigenstates of s°p. This proves 

Design Rule 3.7. Uniaxial spinometers with synoptic or batrachian tunings, but not er- 
godic tunings, asymptotically einselect eigenstates of the measurement generator. 

3.4.4- Convergence bounds for the einselection of eigenstates We now prove a bound on 



the convergence rate of Design Rule 3.7 For QMOR purposes, this bound provides an 
important practical assurance that an ensemble of uniaxially observed spins never becomes 
trapped in a "dead zone" of state space. 

To prove the convergence bound, we notice that in the continuum limit 9 <^ 1 the 



increment (86) can be regarded as a differential equation in simulation time t = nSt. For 
synoptic tuning the inequality £'[A^(s°p)] > £'[A„(s°p)]^ then allows us to derive — by 
integration of the continuum-limit equation — the power-law inequality 



i?[A„(s°P)] < ii;[Ao(s°P)]/(l + 4n02)_ 



(88) 
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This implies that the large-n variance is 0{n^^). This proof nowhere assumes that the ini- 
tial ensemble is randomly chosen; therefore the above bound applies to all ensembles, even 
those whose initial quantum states are chosen to exhibit the slowest possible einselection. 
We conclude that for synoptic tuning the approach to the zero- variance limit is never patho- 
logically slow. We have not been able to prove a similar bound for batrachian tuning, but 
numerical experiments suggest that both tunings require a time t ~ 5t/9'^ to achieve eins- 
election. Proofs of stronger bounds would be valuable for the design of large-scale QMOR 
simulations. 

3.4-5. Triaxial spinometers We now consider triaxial spinometers, in which three pairs of 



synoptic measurement operators (82) are applied, having as generators the spin operators 
{sx,Sy,Sz}, applied with couplings {6j;,6y,9z}. 

3.4.6. The Block equations for general triaxial spinometers In the general case we take 
^1 / ^2 / ^3- We define Xn = {xn,yn, Zn} = j{ipn\s\ipn) to be the polarization vector at 
the n'th simulation step. This vector is normalized such that \xn\ < 1, with \xn\ = 1 if and 
only if \^lJn) is a coherent state. We further define 6xn = Xn — Xn-i- Taking as before E[. . .] 
to be an ensemble average over simulations, such that the density matrix of the ensemble 
is pn = E[\^lJn){1pn\], and therefore = jtispn, we readily calculate that the Bloch 

equation that describes the average polarization of the ensemble of simulations is 



E[5Xn] 
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r q2 
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E[Xn-l] 

E[yn-i] 

E[Zn-l] 



Since it depends only linearly upon p„, the above expression is invariant under the U- 
transform of the Theorema Dilectum. We are free, therefore, to regard our spinometers as 



being ergodically tuned (80), such that the simulation can be equivalently regarded, not as 



three competing axial measurement processes, but as independent random rotations being 
applied along the x-axis, y-axis, and z-axis. The above Bloch equation therefore has the 
functional form that we expect upon purely classical grounds. 

3.4. 7. The einselection of coherent states Now we confine our attention to balanced triaxial 
spinometers, i.e., those having with 61 = 61 = 61 = 6, such that no one axis dominates the 
measurement process. Numerical simulations suggest that for synoptically tuned measure- 
ment processes, in the absence of entangling Hamiltonian interactions, simulated quantum 
trajectories swiftly converge to coherent state trajectories, regardless of the starting quan- 
tum state. We adopt Zurek's (exceedingly useful) concept of einselection |201] to describe 
this process. We now prove that synoptic spinometric observation processes always in- 
duces einselection by calculating a rigorous lower bound upon the rate at which einselection 
occurs. 

Given an arbitrary state we define a spin covariance matrix A„ to be the following 
3x3 Hermitian matrix (of c-numbers): 



(An)fc/ = {ipnlSkSlllpn) " (V"™ I Sfc | ^n) (V'n | | ^n) 



(90) 



This matrix covariance is a natural generalization of the scalar variance An{s°^) (83), and 



in particular it satisfies a trace relation that is similar to ( 84 ) 



tr 



= j if {ipn) is a coherent spin state, 
> j otherwise. 



(91) 
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Here a coherent spin state is any spin-j state \x), conventionally labeled by a unit vector x, 
such that = jx (see, e.g., Perelomov |153| eq. 4.3.35]). The algorithm of Fig. [T] 

evolves the mean spin covariance according to 



^=1 je{(+),(-)} 



{i>n\{M^)^SlM^\tlJn){i'n\{M^)^SmM^\i;n) 



(92) 



(V„|(M;)tM;|V„) 

For compactness we define the A-increment 5 An = A„+i — A„. Then by a series expansion 



of (92) similar to that which led from (85) to (86) — but with more indices — we find that 



for ergodic, synoptic, and batrachian tuning the mean increment is 

{0 ergodic tuning, 

-46'^ trS[A„-A*] synoptic tuning , (93) 
(see text) batrachian tuning. 

The "see text" for batrachian tuning indicates that we have found no closed-form expression 
simpler than several dozen terms; numerical experiments show that for this tuning the 
covariance exhibits random jump- type fluctuations that seemingly have no simple limiting 
behavior. In contrast, synoptic tuning's increment has a strikingly simple analytic form, 
which was guessed as an ansatz and subsequently verified by machine algebra. 

Proceeding as in the proof of Theorem 1, and temporarily omitting the subscript n for 



compactness, we now prove that for A computed from |^) by (90 ), the scalar quantity tr A-A* 
is non-negative for all \^), and vanishes if and only if is a coherent state. We remark 
that this proof is nontrivial because tr A-A* is by no means a positive-definite quantity for 
general Hermitian A; an example is A = [_'^ "''q] , for which tr A-A* = tr [~g _5] = —2. 

Because A is a Hermitian 3x3 matrix, it can be decomposed uniquely into a real 
symmetric matrix A and a real vector v = /j by Aj^ = Kik+i/2 Y^=i ^iki'Vi- Because 



the increment (93) is a scalar under rotations, without loss of generality we can choose a 
reference frame having basis vectors {x, y, z} such that {v-x, v y, v-z} = {0, 0, z} and z >0. 
In this reference frame, the following decomposition is valid for any Hermitian matrix A 
(i.e., it holds for A an arbitrary symmetric matrix and z an arbitrary real number): 

trA-(A*)=j2(l-z2)/2 + /(i_ ^2)2/4 

+ 2pa + ^Ph + IPc + Ipd + IPe + jPU (94) 

where the residual terms Pa.,Ph^ - • • iPf are 

Pa = Ah + A?3 + Als, Pd = (i'(l - Z^) - 2A33)2 

Pb = (All - A22)', Pe = (All + A22)2 - UU + 1) - (A33 + fz^)f 

Pc = (A33)', pf = j2(l-z2)_A33 (95) 

We now prove that each term in this decomposition is non-negative. The terms pa, Ph, Pc, 
and Pd are non-negative prima facie. The term p^ vanishes for arbitrary in consequence 
of the spin operator identity All -I- A22+A33 -I- j^z^ = {i'\sl + sl + sHip) =j{j + l). That the 
remaining terms are non-negative in general follows from the spin operator inequalities —j < 
{'4>\ss\^p) < j and < (V'lsllV') ^ i^) which together with our reference- frame convention 
imply the inequalities — 1 < z < 1 and A33 < — z'^). 
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Next, we show that the sum of terms (94) vanishes if and only if \ip) is coherent, i.e., 



if and only if lip) = \z). It is a straightforward exercise in spin operator algebra to show 
that \ip) = \z) if and only if all of the following are true: z = 1, A33 = j^. An = A22 and 



A12 = A23 = A13 = 0; it follows that (94) vanishes if and only if 1-0) = \z). By reasoning 
similar to Theorem 1, we conclude: 

Design Rule 3.8. Triaxial spinometers with synoptic tunings asymptotically einselect co- 
herent spin states. 

3.4.8. Convergence bounds for the einselection of coherent states We now exploit the iden- 



tity (94) to prove a bound on the convergence of Design Rule 3.8 Our strategy is simi 



lar to our previous proof of the bound for Design Rule 3.7 Substituting the identity 
— z^) = tr A„ — j in the first two terms of (94), and taking into account the non- 



negativity of the remaining terms pg,, ph, Pc, Pd, and pe, we obtain the following quadratic 
inequality in (tr A„ — j): 

trA^-A^ > l{tTAn-j) + l{iTAn-jf>0. (96) 



As an aside, our starting identity (94) was devised so as to imply a general inequality having 
the above quadratic functional form, in service of the proof that follows, but we have not 
been able to prove that the above coefficients {5,3} are the largest possible. 

Upon taking an ensemble average of the above inequality, followed by substitution in 



(93), followed by a continuum-limit integration, we obtain the following convergence bound 



for the ensemble-averaged trace covariance: 

tr^[A„J-j < (trS[Ao] - j){exj>{2n9^) - 1) + 2 exp(2n02) • (^^) 

It is instructive to restate this bound in terms of simulation time t = n5t. Taking the 
continuum limit A„ — > A(t), noting that the timescale Ti = 5t/9'^ is the conventional Ti 
that appears in the Bloch equations ( |89[ ), defining for compactness of notation the initial 
trace covariance to be kq = tr i?[A(0)] — j, and assuming for the sake of discussion that 
kq S> 1 (i.e, we assume that the initial ensemble is far-from-classical) the functional form 
of the above bound exhibits three asymptotic intervals, whose t-dependence is respectively 
0{l/t), and 0(exp(-2i/Ti)): 



{ivE[K{t)]-j) 



< 



Ko(l-Kot/Ti) for 0<t/ri<l/Ko, 

Ti/t for 1/ko<VT'i<1 (98) 

2exp(-2t/ri) for 1 < t/Ti 



The 0(1) and 0{l/t) behavior is functionally similar to the convergence bound established 
in ( |88[ ) for the eigenvalue variance of Theorem 1, namely, an initial linear decrease, followed 
by an 0{\/t) fall-off. Unique to triaxial spinometry (as far as the authors know) is the final 
exponentially rapid convergence to a coherent state. 

We note that convergence is complete within a time ~ Ti that is independent of both 
the spin quantum number j and the overall dimensionality of the Hilbert space in which 
the spin is embedded. As with Theorem 1, this is a worst-case bound that applies to all 
ensembles, including (for example) exotic ensembles initialized with "Schroedinger's cat" 
states. More particularly, it applies to large-dimensional ensembles in which each of n spins 
in a Hilbert space of overall dimension (2j -|- 1)*^ is synoptically observed. 
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3.4-9. Implications of einselection bounds for quantum simulations We now begin to have 
a quantitative appreciation of the geometric assertion of Fig.|4]^i), that quantum simulations 
can be regarded as theaters in which the trajectory compression of synoptic observation op- 
poses the creation of entanglement by Hamiltonian dynamics, with the balance between 
compression and expansion determining the dimensionality of the QMOR state space re- 
quired for accurate simulation. 

Even stronger convergence bounds than those we have proved would be valuable in 
designing QMOR simulations. Especially useful would be more tunings in which noise is 
realized as an entangled measurement. Physically speaking, an entangled measurement is 
performed by interferometrically splitting a photon along n paths, scattering the photon 
from a different spin along each path, then recombining and measuring the photon by freely 
choosing among any of an exponentially large set of braidings and interferometric couplings 
of the downstream optical fibers. 

The analysis of such noise-equivalent tunings would require mathematical methods con- 
siderably more sophisticated than those we have deployed in this article. A known con- 
sequence of the Holevo-Schumacher- Westmoreland (HSW) theorem (which is the quantum 
analog of the Shannon channel capacity theorem) is that entangled measurements are nec- 
essary to maximize the information capacity of quantum channels |il46| sec. 12.3.2]. 

If we hypothesize that quantum trajectory compression is in some sense proportional 
to information extracted by measurement, then the HSW theorem tells us that entangled 
measures will be more effective for QMOR purposes than the single-spin measures that 
we consider in this article. It is likely, therefore, that the search for more efficient QMOR 
techniques will benefit considerably from continued progress in quantum information theory. 

3.4.10. Positive P -representations of the thermal density matrix Now we focus upon con- 



trol and thermodynamics. For i the thermal axis defined in (61), we modify the synoptic 
spinometer matrices such that 

M(*^) = e-^°^(*^«)fe[cos(0Sfc) +sin(0Sfc)]/^, (99) 
Mf_) = e+*"''(*^^)*[cos(0Sfc) -sin(0Sfc)]/^/2, (100) 

where a is the control gain. We will call this a closed-loop triaxial spinometer with unitary 
feedback, because (as we will see) the unitary operators exp(ibia0t x s) act cumulatively to 
align the spin axis with t. 

Closing the control loop does not alter the coherent einselection because the sole effect 
of a post hoc unitary operator on <t„ is a spatial rotation. Since tr <t„-<t* is a rotational 



scalar, (93) still holds. Thus we have 



Design Rule 3.9. Closed-loop triaxial spinometers with unitary feedback asymptotically 
einselect coherent states. 

The density matrix p of an ensemble of closed-loop triaxial spinometer simulations is de- 
scribed by sequence {pi, p2, ■ ■ .} whose increment is 

3 

6pn = 5^ (Mf^)P„MS + Ml^PnM}^^ - pn) , (101) 

fc=i 

By a straightforward (but not short) calculation we find that 5pn vanishes for pn = p^^ if 
and only if the closed-loop gain a satisfies 

a = -tanhi/3 or 1/a = - tanh i/5 . (102) 
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The following two trigonometric identities hold for either choice, and will be used in Sec- 



tion 4.1 to establish that the choice is immaterial in practical numerical simulations. 

1/a + a = — 2 coth and 1/a — a = — 2 csch (103) 

Defining as usual the dimensionless temperature T = l/f3, we see that an optimal control 
gain a — > ±1 establishes a temperature T — > =pO, while a control gain \a\ ^ 1 establishes a 
finite temperature. We will establish later on that p^^ solves 6pn = uniquely, because the 
Fokker-Planck equation for p has a unique stationary solution (thus the approach of the 
density matrix p to thermodynamic equilibrium never "stalls" or becomes trapped at false 
solutions). These results prove 

Design Rule 3.10. The density matrix of an ensemble of closed-loop triaxial spinometer 
simulations is asymptotically thermal. 



To connect ( 101 ) with the thermodynamic literature, we set t = (0, 0, 1) and expand to order 
9'^. The result is equivalent to a thermal model given by Perelomov (eq. 23.2.1 of [153j ). 



Gardiner gives a similar model (eq. 10.4.2 of [84 

^Pn -- 



In Lindblad form we find 



^7(1/ + 1) (s+s_/) 
- l-fiy{s-s+p 



2s_/9S+ + /9S + S_) 

2s+/js_ + ps-S+) 

^ {S^S^P - 2s-ipS3 + S3S3) , 



(104) 



^S2)/v2 are raising and lowering operators, and 
= -406*2 and = -1/2 - (q + l/a)/4. 



where s+ = (si +is2)/V^ and s_ = (si - 
we have adopted Perelomov's variables 7 

3.4.11. The spin-1/2 thermal equilibrium Bloch equations The special case of spin-1/2 
particles in thermal equilibrium often arises in practice. Setting the polarization axis t = z, 
and allowing independent spinometric couplings {0x,9y,0z} as in (89), we find that the 



" E[SXn] ' 


1 

~ ~2 


ia^9l + el + el)E[xn-i] 


E[6yn] 


{el + aX + 0''.)E[yn-i] 


_ E[6zn] _ 




_ {l + a^){dl + el){E[zn^i]+tanhlp) 



finite-temperature synoptic measurement operators (99) imply the following asymmetric 
Bloch equations (valid for j = 1/2 only): 



(105) 



As expected on thermodynamic grounds, we see that the equilibrium polarization is E[z] = 
— tanh 2/5- These equations are a generalization of the usual Bloch equations, in the sense 
that the relaxation rates along the x-, y-, and z-axes can differ independently. We remark 
that for j > 1/2 the thermal Bloch equations do not have a closed analytic form; that is 
why this more general case is not considered here. 

3.4.12. The spinometric ltd and Fokker-Planck equations Now we focus on Ito and Fokker- 
Planck equations, aiming by our analysis to obtain both the already-validated positive P- 



representation of Design Rule 3.2 and (in the large-j limit) both the linear Design Rules 
3.3-5 and the fundamental quantum limits of Design Rule] 3. 6] 



(106) 



We define a binary data three-vector d„ = [d\, d^, d^) by 



+1 
-1 



for \ll}n+l, 
for IV'n+l! 



Then {di, <i2, . . .} is a binary data record with calibration -^[dn] = -E'is^n]) where 



5s = 26 j 



(107) 
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is the spinometer gain. We define a zero-mean stochastic variable Wn by 



(108) 



such that (to leading order in 9) Wn has the second-order stochastic properties of a discrete 
Wiener increment: 

E[{Wn)k{Wn')k'] = Snn'Skk' ■ 



Then via an identity valid for a coherent state, 

{ipnlskSlllpn) = \j5kl + j{i - \)xkXj + \ij eklmXm , 



the spinometer increments (99-b) are equivalent to an Ito increment 

6Xn = Xn+l - Xn= gla{Xn) + Qs h{Xn)-Wn ■ 

For the drift vector a and diffusion matrix b we find 



(109) 

(110) 

(111) 



a{x) 



b{x) 



4i 



^x [a{2j - l)x-i + (1 + la^^ 



+ ^i[a{2j + l)-la^x-i] 



2j 



[I — X ® X + a{i X — x t I)] . 



(112a) 
(112b) 



Design Rule 3.9 asserts that the Ito increment (111) confines the trajectory of Xn to the 
unit sphere. The mean increment of the m'th radial moment \x\'^ must therefore vanish 
when I a; I = 1. We check this by direct calculation, finding 



5En[\x\'^] oc \m{m - 2)E[xn-h{xn)-b\xn)-Xn\ 

+ mE[\Xn\^{Xn-a{Xn) + \ iTb{Xn)-h\xn))] 



(113) 



which indeed vanishes for the a and b of ( |112a -b). By well-known methods [84^, the Ito 
increment (111) immediately yields a Fokker-Planck equation for the PDF Pj{x). Setting 
z = x t we obtain the stationary state equation 







q(1 + z'^) + 2ja(l - z^) - z{l + a^)]pj{z) 



dz 

+ ~[{l-z^){l-2az + a^)]pj{z), 



which (when properly normalized) has a unique solution 

Pj{x) = {a + l/a- 2z)-'^^-'^ 

X (2i + l)[(a + l/a-2)"2j'-i - (a + 1/a + 2)-2i-i] -^/vr 



(114) 



(115) 
(116) 



in which we see that the symmetry a ^ 1/a is indeed respected. Substituting q -|- 1/a = 
— 2coth^/3 per (103), and adjusting the normalization to match the P-representation con- 



vention (59) yields Design Rule 3.2 
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3.4-13. The standard quantum limits to linear measurement To connect these results to 
Design Rules 3.2-3.6 we first write the Ito equation (111 ) in Langevin form by substituting 

5xn Jl^^'dt' x{t') (117) 

a{xn) ^r//+^k'a(a;(t')), (118) 

hix^yWn -^rg,jl+''dt'b{x{t'))-x^{t'), (119) 

where r = l/5t is the rate at which increments occur, and x^{t) is white noise with cross- 
correlation 

E[x'^{t)xl{t')] = 6kk'S{t - t')/{rgl = 6kk'S{t - t')/i4rfe^) . (120) 

(121a) 

where x^{t) = x{t) + x^{t) is the measured spin axis. 

We see that x{t) is dynamically attracted toward the measured axis x^{t). Even open- 
loop spinometers exhibit this attraction, since for a = we find 



Then (111 ) becomes the integral of the Langevin equation 

X = rg^[a{x) + b{x)ix^ 



X 



~x + ^{I 



„M 



x\a=o = rg;^[-^^x + —{I - x(i$ x)-{x"^ - x)] . (121b) 

We remark that in uniaxial spinometry we saw that a similar einselection-by-attraction 
generates the "collapse" of \ipn) to an eigenstate, as described by Design Rule 3.7 This 



attraction is of course a fundamental tenet of the Hilbert ontology of Section 3.2.7 



We now transform (121b) to the second-order Newtonian equation of an oscillator. To 
do this, we introduce a spring k and frequency ujq by defining the operators 

g°P = (hujo/jk)^/'^ [+sicosujot — S2s'mujot) , (122a) 

p°P = {kh/juJaYl'^ (-sisinwoi - S2CosujQt) . (122b) 

We confine our attention to those coherent states that have z ~ — 1, which with our sign 
conventions means systems having positive inverse temperature /5, negative Hilbert feedback 
gain a, and oppositely directed spin x and polarization t, such that x ■ t ~ —1. For 
these states the canonical commutator [g°P,p°P] = —ihs^/j ~ ih holds in the large-j limit. 
Defining the coherent oscillator coordinate q{t) to be 

q{t) = {jhujo/k)^^'^ (x(t) cos Wot — y(t) sincjo^) , (123) 



we find that (121b) takes the linearized Newtonian form 

mq=-kq + f°, (124a) 

q^ = q + q'', (124b) 

Here the spring A;, mass m = kjuj^^ and coordinate q are to be understood in a generalized 
sense in which the system energy is \m(^ -|- }^k(^ . The measurement noise (/"(t) is given 
from (123) in terms of the spinometer noises x^(t) and y^(t) of (120) by 

(125) 



f{t) = {jhujQ/k)^/'^ {x^{t) cos ujQt - y^{t) sin uot) , 



and we find from (120) and (123) that the measurement noise (7°(t) has a PSD Sqn of 

hiOo/{4krje^) (126) 
'y'H[q'^{t)], where Ti is 



Sgn{uj)\ 



The force noise f^{t) is then determined from (121b) to be /°(t) 
the Hilbert transform and the Hilbert gain 7 is found to be 

7 = AkrjO'^/ujQ 



(127) 
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3.4- M- Multiple expressions of the quantum noise limit It follows from the preceding re- 
sults that the PSD of the spinometric force noise f^{t) can be expressed in multiple equiv- 
alent forms: 



expression 

fe2 



h'/S,n{0j)\ 

4khrj9'^/uJo 



physical interpretation 

force noise tx measurement noise (128a) 

force noise oc 1/ (measurement noise) (128b) 

force noise oc Hilbert backaction gain (128c) 

raw spinometer parameters (128d) 



Each of the above relations has a plausible claim to expressing the "most natural" or 
"most fundamental" relation between measurement noise and force noise . . . despite the 
fact that no two physical interpretations are the same, and even though the interpretations 



given (128a) and (128b) seem contradictory. We further see by Design Rule 3.6 that these 



spinometric relations saturate the Hilbert noise limit {jSqn = h) , the Braginsky-Khalili limit 
{SqnSfn = h?), and the Hefner-Haus- Caves limit (NF = 2); thus in some sense all of these 
fundamental quantum limits are embodied in the above family of spinometric relations. 

Acknowledging the self-consistency of this diversity, and appreciating its mathematical 
origin in the diversity of equivalent noise models that are supported by the Theorema 
Dilectum, helps us appreciate how the quantum noise literature can be so immensely large, 
and support so many different notations, physical arguments, and conclusions, and yet 
maintain its internal consistency. 

In a teaching environment, it is not practical to sustain a dispassionately anarchical 
equality among physical interpretations ( 128aj -d). This article's Hilbert ontology (Section 
3.2.4) designates (128a) to be the fundamental relation, because it embodies the central 



Hilbert tenet that "measurement noise always back-acts upon system dynamics in such a 
way as to bring the state of the system into agreement with the measurement." This choice 
is justified solely because it yields useful guiding principles for efficient quantum simulations. 

3.5. Summary of the design rules 



In summary, we have established by Design Rules 3.7 3.8 and 3.9 the quantum mechanism 
by which synoptic noise processes compress simulated quantum trajectories onto lower- 



dimension GK manifolds (as was promised in Section 1.5.11 ). We have established by Design 



Rule 3.10 that the effects of thermal reservoirs can be modeled as equivalent processes 
of covert measurement and control (as was promised in Section 1.5.12). And we have 



established by Design Rule |3. 6] that the Hefner-Haus- Caves, Braginsky-Khalili, and Hilbert 



quantum noise limits are respected by QMOR simulations (as was promised in Section 3.2.6 ). 

The focus of the remainder of this article is to show, by explicit examples, that these 
design rules are sufficient to "enable the reader to proceed to the design and implementation 
of practical quantum simulations, guided by an appreciation of the geometric and informatic 
principles that are responsible for the overall simulation accuracy, robustness, and efficiency" 
(as was promised in the Introduction). 

4. Examples of quantum simulation 

Now we turn our attention toward applying the preceding results in implementing practical 
quantum simulations. 
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4.1. Calibrating practical simulations 



Our simulations provide data via the binary stream of defined in (106), which is low-pass 



filtered to produce a classical data record. We now work through, in detail, the process 
of computing and calibrating this data stream, and ensuring that it is numerically well- 
conditioned. 

We begin by considering the problem of determining, from physical system parameters, 
the measurement operation parameters {6, a} in (99) and the clock rate r = l/5t. In essence 



this calibration process requires that we invert systems of equations that include the Bloch 
equations (89) and (105), the Langevin equation ( 121b[ ), and the mapping of spinometer 
parameters onto oscillator parameters (refeq: q definition). 

Whenever our simulations include a projective step, we must also ensure that the e- 
parameter of the projected Theorema Dilectum (54) satisfies e ^ 1. Physically speaking. 



imposing the small-e condition ensures that the simulated trajectories evolve by drift and 
diffusion, rather than by "quantum jumps" that may be projectively ill-conditioned. Equally 
importantly, it ensures that they respect the informatic causality that is guaranteed by the 
Theorema Dilectum. 



.1.1. 



Calibrating the Bloch equations A common system to be simulated is a spin j = \ 
in contact with a thermal reservoir. We desire that the three thermal relaxation rates 
along the x-, y-, and z- axes be -[T-^;, Ty, T^} — {l/-^; '^l^y^ and that the equilibrium 

thermal density matrix be p^^ oc exp(— /Js^), such that the equilibrium spin polarization po 
is — tanh \l3 as in ( 105 ) . 



We thus have four physical parameters \Tx,Ty,T z, /?} with which to determine five raw 
spinometer parameters {9x,6y,9z,a,r}. Needing one more physical parameter, we note 
that the e-parameters of the GK-projected Theorema Dilectum (54) is given for j = 1/2 
spinometers from ( 105 ) by = 0^(1 -|-a^)/4, and so we impose as our fifth condition e < 0.1 
for all three spinometers (the cut-off 0.1 yielding in our experience well-converged numerical 
results). The equations below then follow from (101), (103), and ( |105 ): 



— tanh or — l/tanh|/? (freely chosen) 
€-= [r, - sgn(l - a^)(Tx - Ty) cosh /(4r) 
e; = [r, + sgn(l - a2)(r, - Ty) cosh /(4r) 

e-= [r, + r,-r,]/(4r) 

?2 = 4e2/(l + «') for 



a 

el 

4 
,2 



(129a) 
(129b) 
(129c) 
(129d) 
(129e) 



Bloch-parameter calibration proceeds as follows. We first determine the spinometer gain 



a from (129a), and we will see that the choice "gain too big" versus "gain too small" is 



immaterial. The value of the spinometer click-rate r is then set from (129b -d) by requiring 
that mm{€x, ey,€z} ^ 0.1 for the reasons noted above. The values of the three spinometer 



phases {Ox,0y,6z} are then determined from (129e 



We remark upon three features. First, we see that insofar as simulation efficiency and 
numerical conditioning are concerned, the choice between the two options for the feedback 



gain a in (129a) is immaterial, since according to the above construction the spinometer 
click-rate r and the e-parameters are unaffected. Second, for the special case Tz = Ti and 
Tx = Ty = T2 the above results reduce to the usual Bloch equations. Third, the positivity of 



the e-parameters in (129b -d) requires that the Bloch relaxation rates satisfy the inequality 



[Tx - r^yl cosh <Tz < {Tx + Ty) 



(130) 
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The authors suspect that the above Bloch inequahty is tight, in the sense that no spin- 1/2 
Lindblad-form master equation can violate it, but we have not proved this. We remark 
that the above triaxially asymmetric Bloch equations and their associated relaxation rate 
inequality have (to the best of our knowledge) not appeared in the literature before. 

4-1-2. Calibrating test-mass dynamics in practical simulations Now we consider test-masses 
{c-g-, MRFM cantilevers) in contact with thermal reservoirs. Calibration proceeds by a line 
of reasoning similar to the above. We take the test-mass to be described by two physi- 
cal parameters: the (dimensionless) temperature /3 = hcoo/^kBT) and the (dimensionless) 
quality Q of the ring-down wave-form q{t) oc cos(ix'ot) exp{uot / {2Q)) . The four spinometer 
parameters to be determined are {j, 6, a, r}. The spin number j we take to satisfy j ^ 1, 
and we will find that the precise value chosen for j is immaterial. For coherent states with 
z ~ — 1 as discussed in Section 3.4.12, the e-parameter is found to be e ~ i^^(l + «^)- 
Calibration proceeds as follows: 



— tanh or — l/tanh^/3 (freely chosen) 
e- = {uJo/{2Qr)) coth i/3 



a 

2 



(131a) 
(131b) 
(131c) 



We first determine the spinometer gain a from (131a). The value of the spinometer click- 
rate r is then set from (131b) by requiring that e < 0.1 as in the Bloch equation case, 
and the spinometer phase 9 is determined from (131c). We remark again that insofar as 



simulation efficiency is concerned, the choice between the two options for the feedback gain 



a in (131a) is immaterial, since the spinometer click- rate r is unaffected. We see also that 



for fixed quality Q, the simulation rate r is ©(coth^/S), which physically means that hot 
cantilevers are computationally more expensive to simulate than cold ones, as is reasonable. 

4-1-3- Calibrating purely observation processes It can happen that we wish to directly 
observe a spin- 1/2 particle along a single axis, nominally the z-axis, in an observation 
process in which the measured z-axis polarization z™(t) = z{t) +z^{t) has a specified (one- 
sided) noise PSD S'^n. No thermo dynamical feedback is applied. Then by reasoning similar 



to the preceding cases, we find from (120) that the calibration relations are 



l/(2r5. 



(132a) 
(132b) 



and as before, the spinometer click rate r is determined by requiring e < 0.1. 

Similarly, if we wish to simulate the continuous observation of a test-mass coordinate 
q{t), with no thermodynamical feedback, the required calibration equations are given from 



(120) and (123) in terms of the measurement noise PSD Sqn by 



(133a) 
(133b) 



where again the spinometer click rate r is determined by requiring e < 0.1. 



We remark that in all of the above cases the raw binary data stream ( 106 ) of the simu- 



lation must be low-pass filtered in order to obtain a (noisy) data record that will have the 
above statistical properties within the filter passband having. This filtering closely models 
the way that real experimental signals (for example, a continuously measured cantilever 
motion) are displayed upon oscilloscopes. 
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remarks on batrachian unraveling 

• quantum jumps thermalize the spin 

• measured rms polarization is unity 

• quantum rms polarization is unity 



remarks on ergodic unraveling 

random fields thermalize the spin 
measured rms polarization is unity 
quantum rms polarization is classical 



remarks on synoptic unraveling 

• measurement noise thermalizes the spin 

• measured rms polarization is unity 

• quantum rms polarization is subclassical 



Figure 10: Simulation of single electron moment detection by MRFM. 



4-2. Three single-spin MRFM simulations 

With reference to Fig. [TOj we now turn our attention to the simulation of the IBM single- 
spin MRFM experiment |il70j . We will initially present the simplest possible class of 
simulations that reproduce the data of that experiment, postponing a discussion of more 



detailed simulations until Section 4.4 Our goal is to illuminate the central role of the 
Theorema Dilectum in answering the question that was raised in Section |1.1[ "How does 
the Stern-Gerlach experiment work?" 

All three columns of Fig. 10 show a simulated thirteen-hour experiment (the length of 
the IBM experiment). The time-spacing 6t = 1/r between spinometric clicks is set to 7.1 ms; 
thus approximately 6.6 x 10^ time-steps were simulated. In each column, the experimental 
data are simulated as arising from three competing spinometric processes. Spin relaxation 
was simulated by x-axis a nd y -axis spinometers having 6x = 6y = 0.093. The consequent 
spin relaxation time from (105) is Tz = 2/(r(0^ -|- e^)) ~ 0.76 s as was observed in the IBM 
experiment. Measurement effects were simulated by a z-axis spinometer having 6z = 0.026, 
and the consequent measurement noise PSD from (120) is Sz^ = 2/{r0l), which numerically 
corresponds to a noise level of 11.5 Bohr magnetons of noise in one root-Hertz of bandwidth, 
as was observed in the IBM experiment. 

For visualization purposes only, all time-domain data streams shown in Fig. 10 were 
low-pass filtered with a time constant r = = 0.76 s. 

The following discussion is insensitive to the above experimental details, and applies to 
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all experiments of Stern-Gerlach type in which the signal-to-noise ratio (SNR) is low, such 
that continuous monitoring over extended periods of time is required to observe the effect. 

In the next three sections, we simulate the spin relaxation of the IBM single-spin MRFM 
experiment by three different unravelings: batrachian, ergodic, and synoptic. We will see 
that the three unravelings lead to three very different classes of quantum trajectories, and 
hence, three very different answers to the question "How does the Stern-Gerlach experi- 
ment work?" Nonetheless, as guaranteed by the Theorema Dilectum, we will find that the 
simulated experimental data are identical for all three unravelings. We now work through 
the mathematical and physical details of how this comes about. 

4.2.1. A batrachian single-spin unraveling The left-hand column of Fig. 



10 



shows a sim- 
ulation in which thermal noise is unravelled as a batrachian process, whose measurement 
operations are given algebraically in (81) and which are depicted in hardware-equivalent 



form in Figure [9[b). This is by far the easiest simulation to analyze in closed form: the spin 
polarization jumps randomly between ±1, driven by the batrachian jumps of the thermal 
reservoir, while being continuously measured by the (noisy) cantilever. 

The simulated data stream is therefore a random telegraph signal with added white 
noise, such that the mean-square quantum spin polarization inferred from the data is unity. 
We conclude that from the batrachian point of view, the Stern-Gerlach effect (meaning, 
that the mean-square spin polarization is measured to be unity) comes about because noise 
is a quantized jump process, such that the mean-square spin polarization always is unity. 



4.2.2. An ergodic single-spin unravelling The middle column of Fig. 10 ^b) shows a sim- 
ulation in which thermal noise is unravelled as an ergodic process, whose measurement 
operations are given algebraically in (80) and which are physically depicted in Figure [9|^a) . 
Physically speaking, the spin polarization is driven by random magnetic fields, such that 
the mean-square quantum polarization is 1/3. 

Now a subtle effect comes into play. The 2;-axis measurement process back-acts upon the 
spin state, such that whenever an "up" fluctuation in the data is observed, the spin state is 
"dragged" toward a positive polarization. This effect is evident in the simulated data. The 
consequence of state-dragging back-action is that the measured mean-square polarization 
is larger than the mean-square polarization of the underlying quantum state. 

It would be quite a complicated task to calculate the resulting data statistics from {e.g.) 
the appropriate Ito, Langevin, and Fokker-Planck equations. Fortunately, the Theorema 
Dilectum does this mathematical work for us: the data statistics are guaranteed to be 
exactly the same random telegraph statistics as in the Batrachian case. 

We conclude that from the ergodic point of view, the Stern-Gerlach effect (meaning, that 
the mean-square spin polarization is measured to be unity) comes about because measure- 
ment is a Hilbert process (meaning, it accords with the state-dragging Hilbert back-action 
ontology of Section 3.2.4). 



4.2.3. A synoptic single-spin unravelling The right-hand column of Fig. 10 'c) shows a 
simulation in which thermal noise is unravelled as an synoptic process, whose measurement 
operations are given algebraically in (82) and which are physically depicted in Figure [oj^c). 

In synoptic unravelling, all processes are measurement processes, and each process seeks 
to align the spin polarization along its own axis. In our simulation, the x-axis and y-axis 
measurement processes are considerably stronger than the z-axis process. In consequence 
of the Hilbert state-dragging effect, the spin polarization now points predominantly in the 
equatorial direction, such that the mean-square quantum polarization is only ~ 0.05. 

Again it would be quite a complicated task to calculate the resulting data statistic from 
Ito equations, etc., and again the Theorema Dilectum does this mathematical work for us: 
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as in the preceding two cases, the data statistics are random telegraph statistics with added 
white noise. We conclude that from the synoptic point of view, the Stern-Gerlach effect 
is not associated with "wave function collapse," but rather comes about (as in the ergodic 
case) because measurement is a Hilbert process. 

4-3. So how does the Stern-Gerlach effect really work? 

We are now in a position to answer more completely the question "How does the Stern- 
Gerlach effect really work?" We answer as follows: "Nothing definite can be said about the 
internal state of noisy systems, either at the classical or at the quantum level. It is best to 
pick an ontology that facilitates rapid calculations and suggests interesting mathematics. 
For purposes of large-scale quantum simulation, a particularly useful ontology is one in 
which all noise processes are conceived as equivalent covert measurement processes. In 
this ontology, the Stern-Gerlach effect works because competing measurement processes 
exert a Hilbert back-action mechanism that 'drags' quantum states into agreement with 
measurement. In consequence of these competing Hilbert measurements, experimental data 
having the statistics of random telegraph signals are obtained even when no quantum jumps 
are present." 



Of course, we saw in the three simulations of Fig. lOi-c that other explanations are 



perfectly reasonable, and that is why in Section 1.1.1 we embraced Peter Shor's maxim: 
"Interpretations of quantum mechanics, unlike Gods, are not jealous, and thus it is safe to 
believe in more than one at the same time," to which we now append the caveat "provided 
that all interpretations respect the fundamental mathematical and physical invariance of 
the Theorema Dilectum." 

4-4- Was the IBM cantilever a macroscopic quantum object? 

The least realistic element of the proceeding simulations is the modeling of the cantilever 
as a single z-axis spinometer having quantum number j = 1/2. A more realistic model 
would have treated the cantilever as a large-j quantum object subject both to thermal 
noise processes and to experimental measurement processes. However, we can appeal to 
the Theorema Dilectum to show that these refinements would not change the simulated 
data at all. The reason is that both the cantilever thermal reservoir and the experimen- 
tal (interferometric) cantilever measurement process can be modeled as synoptic processes 
that compress the cantilever's quantum state to a coherent state. Then modeling the spin 
relaxation as a batrachian process, the output of the resulting (effectively semi-classical) 
batrachian simulation will be precisely the random telegraph signal that was obtained in 



the simpler batrachian simulation of Section 4.2.1 above 



It follows by the Theorema Dilectum that all quantum simulations of the cantilever, even 
elaborate large-j simulations in which a non-coherent quantum cantilever state is entangled 
with the quantum spin state, will simulate the same random telegraph data statistics as 
the simpler simulations already given, and in particular, will yield an observed mean-square 
polarization of unity. 

This leads to an interesting question: what was the "real" quantum state of the IBM 
cantilever? We have seen that this question has a well-posed answer only insofar as there 
is agreement upon the "real" noise and observation processes acting upon the cantilever, 
such that the tuning ambiguity of the Theorema Dilectum does not come into play. 

If we stipulate that the "real" cantilever thermal noise and the "real" spin relaxation are 
due to ergodic physical processes, then the IBM experiment can only be "really" described 
in terms of a spin-state that is quantum-entangled with the cantilever state, in which the 
observed mean-square polarization of unity is "really" due to the state-dragging Hilbert 
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back-action associated with the cantilever measurement process. In other words, the IBM 
experiment "really" observed the cantilever to be a macroscopic quantum object. 

As quantum objects go, the IBM cantilever was exceptionally large |170] : its resonant fre- 
quency was LJo/(27r) = 5.5 kHz, its spring constant was k = 0.011 mN/m, and its motional 
mass was m = k/ujQ = 9.1 pg. The preceding paragraphs are an argument for regarding 
this cantilever to be among the stiffest, slowest, most massive dynamical systems whose 
quantum nature has been experimentally confirmed. Such measurements are significant 
from a fundamental physics point of view, in probing the limits of quantum descriptions of 
macroscopic objects, as reviewed by Leggett |125l 11261 11271 1128j . 

Any line of reasoning that is as brief as the preceding one, about a subject that is as 
subtle as macroscopic quantum mechanics, is sure to have loopholes in it. A major loophole 
is our modeling of decoherent noise as a Markovian process. As reviewed by Leggett et al. 
|129] . spin decoherence in real experiments is (of course) due to non-Markovian quantum- 
entangling interactions. We now turn our attention to the algorithmic and numerical chal- 
lenges of simulating such systems. 

4-5. The fidelity of projective QMOR in spin-dust simulations 

As test cases, we computed what we will call spin-dust simulations. Spin-dusts are quantum 
systems that are deliberately constructed so as to have no symmetries or spatial ordering. 
Their sole purpose is to provide a well-defined test-bed for numerical and analytic studies 
of the fidelity of projective quantum model order reduction. 

Spin-dusts couple pairs of spin-1/2 particles {j, k} via a dipole-dipole interaction Hamil- 
tonian Hjk that is given by 

\ Sj ■[! -2>njk®rijk]- Sk for j / /c 

[Sj-Ujk t0TJ = k 

The unit vectors rijk are chosen randomly and independently for each {j, k}, and we note 
that self-coupling is allowed. Physically we can think of spin-dusts as broadly analogous 
to — but less structured than — systems such as the interacting spins in a protein molecule. 

In our simulations each spin is randomly coupled to four other spins, in addition to its 
self-interaction. Then it is easy to show that triJ = and ti H"^ / dim'H = nspm, which 
is to say, the per-spin energy of our spin-dusts has zero mean and unit variance. The 
time-scale of the spin dynamics of the system is therefore unity. We further stipulate that 
each spin is subject to a triaxial spinometric observation process having relaxation time 
Tx = Ty = Tz = 10. Thus the time-scale of decoherent observation is ten-fold longer than 
the dynamical time-scale. 

Simulations were computed with a time-step 6t = 0.1 and spinometric couplings 9^ = 
^j/ = = 0.1, using the sparse matrix routines of Mathematica. The numerical result was 
an "exact" (meaning, full Hilbert space) quantum trajectory \ipQ{t)). These trajectories 
were then projected on GK manifolds of various order and rank by the numerical methods 
of Section |2.12 The main focus of our numerical investigations was the fidelity of the 



projected states lipjcit)) relative to the exact states \ipo{t))- 



11' 



4-5.1. The fidelity of quantum state projection onto GK manifolds With reference to Fig. 
simulations were conducted with numbers of spins n G 1, 18, having random dipole coupling 
links as depicted. The median quantum fidelity was then computed, as a function of n, for 
GK rank r G {1,2,5,10,20,30} (see Fig. [sjfor the definition of GK rank). Both synoptic 
and ergodic unravelings were simulated. Typically |^o)(^) was projected at thirty different 
time-points along each simulated trajectory, always at times t > 100 to ensure that memory 
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Figure 11: The dependence of QMOR fidelity upon GK order and rank. 



of the (randomly chosen) initial state was lost. We remark that numbers of spins n > 18 
could not feasibly be simulated on our modest computer (an Apple G5). 

The quantum fidelity of a projected state IV'a:) was defined to be [146,, Section 9.2.2] 



/ = K^^IV'o)l/[(V'^l^^)(^olV'o)] 



1/2 



(135) 



As shown in Fig. [TT|b), for ergodic unravellings large-n quantum fidelity fell-off exponen- 
tially, while for synoptic unravelings large-n fidelity remained high. 

No mathematical explanation for the observed exponential fall-off in ergodic unravelling 
fidelity is known. The asymptotic large-n behavior of the synoptic fidelity also is unknown. 
In particular, for systems of hundreds or thousands os spins, would the empirical rule-of- 
thumb "GK rank fifty yields high fidelity for spin-dust systems" still hold true? These are 
important topics for further investigation. 

The achieved high-fidelity algorithmic compression was large: an 18-spin exact quantum 
state \'ipo{t)) is described by 2^^ independent complex numbers, while an order-18 rank 30 
GK state \ipic{t)) — as seen at lower right in Fig. [IT|b) — is described by 30 x (18 + 1) = 570 
independent complex variables. The dimensional reduction is therefore 460-to-l. 

4-5.2. The fidelity of spin polarization in projective QMOR We next turned our attention 
to measures of local quantum fidelity, as depicted in Fig. 12 All simulated trajectories in 
this figure were for n = 15 spin-dust. The first such measure we consider are the direction 
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Figure 12: Measures of projective fidelity for n = 15 spin-dust. 



cosines {ip\s ■ rh\ip) / j for randomly spins, randomly chosen trajectory points, and randomly 
chosen unit vectors m. One hundred randomly chosen data points are shown. We observe 
that the rank-one GK manifold does an excellent job of representing the spin direction 
cosines, which can be regarded as (essentially) classical quantities. 

4-5.3. The fidelity of operator covariance in projective QMOR As a measure of pair- wise 



quantum correlation, we examined the spin operator covariance. With reference to (90), 
this quantity is given by 

Sj-fc = 4Aki = {4^Wjak\tp) - (V'|cTj|V)(V'kfe|V'>- (136) 

which vanishes for rank-1 (product states). 

The second row of Fig. [12] plots T,ki for one hundred randomly chosen trajectory points, 
and randomly chosen spins, having randomly chosen indices j and k. We observe that 
GK ranks in the range 20-50 are necessary for projection to preserve pair-wise quantum 
correlation with good accuracy. 

4.5.4- The fidelity of quantum concurrence in projective QMOR As a measure of pairwise 
quantum entanglement, we examined Woofers' quantum concurrence |199j . The concur- 
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rence is computed as follows. Let /Oab be the reduced density matrix associated with spins A 
and B. Let Aj be eigenvalues of the non-Hermitian matrix /OabPab, in decreasing order, where 
Pab = {o'y cry)p*g{ay (8> CTy)- Then the concurrence c is defined to be 



Ai - VA2 - VAs - VA4 (137) 

It can be shown that the concurrence vanishes for product states, and that spins are pairwise 
entangled if and only if c > 0. 



In the third row of Fig. 12 we observe that coupled spin-pairs are far more likely to 
be quantum-entangled than non-coupled spin-pairs, as expected on physical grounds. We 
further observe that GK ranks in the range 20-50 are necessary for projection to preserve 
concurrence with good accuracy. 

4-5.5. The fidelity of mutual information in projective QMOR As a measure of pair-wise 
quantum information, we examined von Neumann's mutual information [ 146] Section 11.3], 
which is computed as follows. For a general density matrix p we define the von Neummann 
entropy S{p) = — tr /olog2 p. Then for spins A and B the mutual information is given by 

S{p,) + S{ps) - S{p,s) (138) 

It is known that the mutual information vanishes for product states, and that this quantity 
is otherwise positive in general. 

In the fourth row of Fig. [T2| we observe that coupled spin-pairs share more mutual 
quantum information than non-coupled pairs, as expected on physical grounds. We further 
observe that GK ranks in the range 20-50 are necessary for projection to preserve mutual 
quantum information with good accuracy. 

As a quantitative summary of this observation, the 15-spin simulations of Fig. |12| predict 
15 single-spin density matrices p^ and 105 pairwise reduced density matrices Pab (in addition 
to higher-order correlations). Each of the single-spin density matrices has 3 (real) degrees of 
freedom, and each pairwise density matrix introduces 9 more (real) independent degrees of 
freedom, for a total of 990 independent degrees of freedom associated with the one-spin and 
two-spin reduced density matrices. In comparison, the rank-50 GK manifold onto which 
the quantum states are projected is described by Kahlerian coordinates having (it can be 
shown) 1600 locally independent coordinates. Using 1600 state-space coordinates to encode 
990 physical degrees of freedom represents a level of MOR fidelity that (obviously) cannot 
be improved by more than another factor of two or so. The mathematical origin of this 
empirical algorithmic efficiency is not known. 

4.6. Quantum state reconstruction from sparse random projections 

We will conclude our survey of spin-dust simulations with some concrete calculations that 
are motivated by recent advances in the theory of compressive sampling (CS) and sparse 
reconstruction. It will become apparent that synoptic simulations of quantum trajectories 
mesh very naturally with CS methods and ideas. To the best of our knowledge, this is the 
first description of CS methods applied to quantum state-spaces. 

Our analysis will mainly draw upon the ideas and methods of Donoho [61] and of Candes 
and Tao [37], and our discussion will assume a basic familiarity with these and similar CS 
articles [SO] EU |32l [311 EB|, especially a recent series of articles and commentaries on the 
Dantzig selector [TTl [291 ESI IMl IZHl 11371 I167j . Our analysis can alternatively be viewed 
as an extension to the quantum domain of the approach of Baraniuk, Hegde, and Wakin 
|111 Wl\ I194j to manifold learning |188| from sparse random projections. 

Our objectives in this section are: 
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• establish that synoptically simulated wave functions tpQ are compressible objects in the 
sense of Candes and Tao ^37j , 

• establish that high-fidelity quantum state reconstruction from sparse random projec- 
tions is algorithmically tractable, 

• describe how nonlinear GK projection can be described as an embedding within a 
larger linear state-space of a convex optimization problem, and thereby 

• specify algorithms for optimization over quantum states in terms of the Dantzig se- 
lector (a linear convex optimization algorithm) of Candes and Tao [37] . 

At the time of writing, the general field of compressive sensing, sampling and simulation 
is evolving rapidly — "Nowadays, novel exciting results seem to come out at a furious pace, 
and this testifies to the vitality and intensity of the field" [38^ — and our overall goal is 
to provide mathematical recipes by which researchers in quantum sensing, sampling and 
simulation can participate in this enterprise. 

4-6.1. Establishing that quantum states are compressible objects To establish that i/^q is 
compressible, it suffices to solve the following sparse reconstruction problem. We begin by 
specifying what Donoho and Stodden [62] call the model matrix and what Candes and Tao 
[37j call the design matrix to be an n x p matrix X. The projected state cf)^ = Xi/jq is given, 
and our reconstruction task is to estimate the ipQ (the "model") from the cf)^ (the "sensor 
data"). In general n < p and we particularly focus upon the case n <^ p. 

We initialize the elements of the design matrix X with i.i.d. zero-mean unit-norm (com- 
plex) Gaussian random variables. Then the rows and columns of X are approximately pair- 
wise orthogonal, such that X satisfies the approximate orthogonality relation XX^ ~ pi 
and therefore satisfies the approximate projective relation (X^X)^ ~ pX^X. As a re- 
mark, if we adjust X to make these orthogonality and projective relations exact instead 
of approximate — for example by setting all n nonzero singular values of X to unity — our 
sparse reconstructions are qualitatively unaltered. 

In CS language, we have specified random design matrices X that satisfy the uniform 
uncertainty principle (UUP) |33| IMl |37] . meaning (loosely) that the columns of X are 
approximately pairwise orthogonal. See [37] for a definition of UUP design matrices that is 
more rigorous and general. 

From a geometric point of view, this means we can regard X'^X — which will turn out 
to be the mathematical object of interest — as a projection operator from our (large) p- 
dimensional quantum state-space onto a (much smaller) n-dimensional subspace. 



We have already seen in Sections 2.12| and |4.5| that the following minimization problem 



can be tractably solved by steepest-descent methods: 

min||V'o-V'.(c)||2. (139) 

c " ^ 

where we have adopted the CS literature's practice of specifying the h norm explicitly. Here 
■j/'k(c) is a vector of multilinear gabion-Kahler (GK) polynomials as defined in Section [2] 
and depicted in Fig. |3j Inspired by the CS literature, we investigate the following CS 
generalization of ([139]): 

mm||X(V'Jc) -^o) =min||<^o-XV'K(c)||^^ (140) 

Now we are minimizing not on the full Hilbert space, but on the n-dimensional subspace 
projected onto by X. We recognize the right-hand expression as a nonlinear Kahlerian 
generalization of a standard minimization problem (it is discussed e.g. by Donoho and 
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Stodden [62l eq. 3] and by Candes and Tao [371 eq. 1.15]). To make this parallelism more 
readily apparent, we can write the above minimization problem in the form 



mm\\y-Xp\\l s.t. /3 = ^,(c) (141) 



for some choice of c, where we have substituted (f)Q ^ y and introduced (3 as an auxiliary 
variable. Comparing the above to the well-known LASSO minimization problem |62| [79] 

min||y -X/3||^^ s.t. < t (142) 



/3 



for some t, we see that the sole change is that the LASSO problem's li sparsity constraint 
—t has been replaced with the GK represent ability constraint /? = i/'^(c). We remark 
upon the parallelism that both constraints are highly nonlinear in f3. 

But this parallelism in itself does not give us much reason to expect that the minimiza- 
tion (141 ) is tractable, since we saw in Section |2] that the space of feasible solutions il^f^{c) 



is (floridly) nonconvex. Consequently, unless some "GK magic" of comparable algorithmic 
power to the well-known "^i magic" of CS theory [35^ should come our rescue, there seems 



to be little prospect of computing the minimum (141 ) in practice. 

Persisting nonetheless, we compute successive approximations {ci, C2, . . . , Cj} by a pro- 
jective generalization of the same steepest-descent method that produced the results of 



Figs. 11 and 12 Specifically, we expand the GK coordinates via Cj+i = Ci + 6ci and iterate 



the resulting linearized equations in 6ci 

6ci = -{A^X^XAf A^X^X{iP^{ci) - i^o). (143) 

Empirically, good minima are obtained from ©(dim/C) iterations of this equation from 
randomly-chosen starting-points. This benign behavior is surprising, given that our objec- 



tive function (140) is a polynomial in ©(dim/C) variables having ©((dim'H)^) independent 
terms, because generically speaking, finding minima of large polynomials is computationally 
infeasible. 



According to the geometric analysis of Section 2.5 the existence of feasibly computed 



minima is explained by the rule structure of GK state-space, which ensures that almost all 



state-space points at which the increment ( 143 ) vanishes are saddle points rather than local 



minima, in consequence of the nonpositive directed sectional curvature that is guaranteed 



by Theorem 2.1 



We now discuss GK geometry from the alternative viewpoint of CS theory, further 
developing the idea that GK rule structure provides the underlying geometric reason why 
CS "works" on GK state-spaces. 

4-6.2. Randomly projected GK manifolds are GK manifolds With reference to the algo- 



rithm of Fig. [2] we immediately identify (A^X^XA)"^ in (143) as the Kahlerian metric of a 
GK manifold having an algebraic Kahler potential (see (|17[)) that is simply 

k{c,c) = ltjj{c)X^X'iP{c). (144) 

Since X is constant, we see that the projected Kahler potential is a biholomorphic poly- 
nomial in the same variables and of the same order as the original Kahler potential. It 
follows that a projected GK manifold is itself a GK manifold, and in particular the GK rule 
structure is (of course) preserved under projection, and this means that all of the sectional 
curvature theorems of Section [2] apply immediately to QMOR-CS on GK state-spaces. 
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This GK inheritance property is mathematicahy reminiscent of the inheritance proper- 
ties of convex sets and convex functions, and it suggests that a calculus of GK polynomials 
and manifolds might be developed along lines broadly similar in both logical structure and 
practical motivation to the calculus of convex sets and convex functions that is presented 



in the standard textbooks of CS [21] (we discuss this further in Section 4.6.9). 

4-6.3. Donoho- Stoddard breakdown at the Candes-Tao bound Putting these ideas to nu- 
merical test, using the same spin-dust model as in previous sections, we find that random 
compressive sampling does allow high-fidelity quantum state reconstruction, provided that 



the state trajectory to be reconstructed has been synoptically unraveled (Fig. 13). 

These numerical results vividly illustrate what Donoho and Stodden [62 have called 
"the breakdown point of model selection" and we note that Candes and Tao have described 
similar breakdown effects in the context of error-correcting codes [33]. Surprisingly, it does 
not appear to have been recognized that a similar breakdown occurs in quantum modeling 
whenever too many wave function coefficients are reconstructed from too few projections. 



For discussion, we direct our attention to the rank-30 block of Fig. 13 ^f), where (as 
labeled) we reconstruct the p = dimTC = 2"=p'" = 2048 (complex) components of tpQ 
from n random projections onto a GK manifold having S = dim/C = (GK rank) x (1 -|- 
log2(dim'H)) = 30 X 12 = 360 (complex) dimensions. All six blocks (a-f) of the figure are 
similarly labeled with p and S. 

Given our example with p = 2048 and S = 360, what does CS theory predict for the 
minimum number n of random projections required for accurate reconstruction? According 
to Candes and Tao 1371 



With overwhelming probability, the condition [for sparse reconstruction] holds for S = 
0(^n/log{p/n)). In other words, this setup only requires 0[log{p/n)) observations per 
nonzero parameter value; for example, when n is a nonnegligible fraction of p, one only 
needs a handful of observations per nonzero coefficient. In practice, this number is 
quite small, as few as 5 or 6 observations per unknown generally suffice (over a large 
range of the ratio p/n). 

We naively adapt the above Candes-Tao big-O sampling bound to the case at hand by 
recalling that S = n/\og{p/n) implies n ~ S\og{p/ S) for p^ n |52]. We therefore expect 
to observe Donoho-Stodden breakdown at a (complex) projective dimension nsh{p, S) (which 
we will call the sampling bound) that to leading order in S/p is 

n,y,{p,S) 2. Slogip/S)\g^^.^j^ (145) 

p = dim Tl 

Here we recall that dimTC is the (complex) dimension of the Hilbert space within which 
the efflorescent GK state-space manifold of (complex) dimension dim/C is embedded (see 
Sections 2.6.4 and 2.11 and also for discussion of how to calculate dim/C). 



The above sparsity bound accords remarkably well with the numerical results of Fig. 13 
This empirical agreement suggests that QMOR and CS may be intimately related, but on 
the other hand, there are the following countervailing reasons to regard the agreement as 
being possibly fortuitous: 

1. the Candes-Tao bound applies to state-spaces that are globally linear, whereas we 
are minimizing on a GK state-space that is only locally linear, and 



2. the onset of Donoho-Stodden breakdown in Fig. 13 is (experimentally) accompanied 



by the onset of multiple local minima of (140), which are not present in the convex 
objective function of Candes and Tao, and 
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Figure 13: Quantum state reconstruction from sparse random projections. A (typical) state 
from an 11-spin trajectory was reconstructed from sparse projections onto random sub- 
spaces (horizontal axis), and the resulting quantum fidelity was evaluated (vertical axis). 



Each point represents a single minimization of (143), by iteration of (143) with a conjugate 
gradient correction, from a random starting point chosen independently for each mini- 
mization. Convergence to "false" local minima was sporadically encountered for low-rank 
GK projections (graphs a-b, GK ranks 1 and 2) but not for higher-rank GK projections 
(graphs c-f, GK ranks 5, 10, 20, and 30). The onset of Donoho-Stodden breakdown was 
observed to occur near the Candes-Tao bound ( 145 ) — plotted as a dotted vertical line — for 



all GK ranks tested. The ergodic spin-dust simulation yielded states whose reconstruction 
properties were indistinguishable from random states, as expected. 
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3. high- accuracy numerical agreement with a "big-O" estimate is fortuitous; the agree- 



ment seen in Fig. 13 is better than we have a reason to expect. 

So although the Candes-Tao bound seems empirically to be the right answer to "when 
does Donoho-Stodden breakdown occur in quantum model order reduction?" the numerical 
calculations do not explain why it is the right answer. 

We now present some partial results — which however are rigorous and deterministic 
insofar as they go — that begin to provide a nontrivial explanation of why the Candes-Tao 
bound applies in the sparse reconstruction of quantum states. The basic idea is to embed 



the nonlinear minimization ( 139 ) within a larger-dimension problem that is formally convex. 
We will show that this larger-dimension optimization problem can be written explicitly as 
a Dantzig selection. 

The main mathematical tool that we will need to develop is sampling matrices X whose 
(small) row dimension is n = dim H, and whose (large) column dimension p is a power of 
dim H. These matrices are too large to be evaluated explicitly — they are what Cai and Lv 
call "ultrahigh-dimensional" [29] . A novel aspect of our analysis is that we construct these 
matrices deterministically, such that their analytic form allows the efficient evaluation of 
matrix products. 

4-6.4- Wedge products are Hamming metrics on GK manifolds Let us consider how the 



efflorescent GK geometry that we described in Sections|1.5.8 and 2.11 can be made the basis 



of a deterministic algorithm for constructing good sampling matrices. 

Our basic approach is to construct a deterministic lattice of points on GK manifolds, 
together with a labeling of the lattice for which the Hamming distance between two labels 
is a monotonic function solely of the wedge product between that pair of points, such 
that the larger the Hamming distance between two points, the closer they approach to 
mutual orthogonality. The problem of constructing good sampling matrices then becomes 
equivalent to the problem of constructing good error correction codes. 

We now construct the desired GK lattice. We consider sampling matrices X whose 
columns are not random vectors, but rather are constrained to satisfy ■0 = ipi^{c) for some 
GK state-space ■j/7^(c). We wish these vectors to be approximately orthogonal. To construct 
these vectors (and simultaneously assign each vector a unique code- word), we specify an 
alphabet of four characters {a, b, d, e}, we identify the four characters with the four vertices 
of a tetrahedron having unit vectors {ria, n^, n^;, rie}, and we identify the unit vectors with 
the four spin-j coherent states {{na), \nb), \nd), I'^e)}, such that for s = {sx,Sy,Sz} the 
usual spin operators, the tetrahedron vertices are ha = {ha\s\ha) /j , etc. Soon it will 
become apparent that the vertices of any polytope, not only a tetrahedron, suffice for this 
construction, and that the vertices of Platonic solids are a particularly good choice. 



We recall from our study of GK geometry that a wedge product (12) can be associated 
to each letter-pair {a, 6} as follows 

\a A = {ha\na){hb\hb) - {ha\hb) {hb\na) (146) 



From Wigner's identity (56) we have \{ha\hb)\'^ = \Djj{0,dab,0)\'^ = cos{9ab/2)^^ where 
cos{6ab) = ha • hh, so the spm-j wedge product is easily evaluated in closed form as 

|aA6|2 = l-cos(W2)^^' (147) 
which for our tetrahedral alphabet is simply 



for a = b 
for a ^ b 



l«A5|,^ = ^ . (148) 
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Here and henceforth we have added a subscript j to ah wedge products for which an analytic 
form is given that depends explicitly on the total spin j. 

Now we specify a dictionary to be a set of n-character words {w^} with each word 
associated with an ordered set of tetrahedral characters w'^ = {c^, C2, • • • , c^}. We further 
associate with each word a petal-vector \w^) 

= |n^,)0|n^,)...®|n,,) (149) 

Given two words {w^,w''}, their mutual Hamming distance h{w'^,w^) is defined to be the 
number of symbols that differ between w'^ and w^. We also can associate with any two 
petal-vectors their mutual wedge product A defined by 

Au-Y = {w'\w'){w''\w'') - {w'\w''){w''\w') (150) 

The Hamming distance h{w'^,w^) is a Hamming metric on our codeword dictionary, and it 
is easy to show that this code metric is related to the petal- vector wedge product \w^ Aw''\ 
by the simple expressions 

^spin 

A tD^^P = 1 — \{h^i^\h^k^)\'^ in general, from 

m=l 

Aw''\'j = I - 9--5'''^J tetrahedral dictionary, from 

or equivalently for a tetrahedral petal-vector dictionary 

hj{w\w'') = -logg{l - \w' A w^\'j)/j (153) 

The main result of this section is the above monotonic functional relation between a Ham- 
ming distance and a wedge product. We are not aware of previous CS work establishing 
such a relation. 

The practical consequence is that given a dictionary of n-character words {w^} having 
mutually large Hamming distances {i.e., a good error correcting code), this construction de- 
terministically specifies a set of nearly-orthogonal petal-vectors (i.e., good vectors for 

sparse random sampling) . Conversely, the general problem of constructing a deterministic 
set of nearly-orthogonal petal-vectors is seen to be precisely as difficult as deterministically 
constructing a good error correcting code. 

4.6.5. The n and p dimensions of deterministic sampling matrices The column dimension 
p of the petal- vector sampling matrices thus constructed is given in Table [T] for Hamming 
distances 1-4 as a function of the row dimension n, the number of polytope vertices n^cv, 
and the dimensionality of the polytope space dim Vq {e.g., for our tetrahedral construction 
nver = 4 and dim Vq = 2). We see that for fixed row dimension n, larger Hamming distances 
are associated with smaller column dimension p, as is intuitively reasonable: the more strin- 
gent the pairwise orthogonality constraint, the smaller the maximal dictionary of sampling 
vectors that meet this constraint. 

The construction has a further dimensional constraint as follows: it is straightforward 
for values of n that are powers of two (because the tetrahedral construction can be used) , 
more complicated when n has factors other than 2 (because larger- dimension polytope 
vertices must be specified), and infeasible when n is large and prime. These constraints are 
reminiscent of similar constraints that act upon the fast Fourier transform, and arise for 
basically the same number-theoretic reason. 



149-150) 



fl511 



147-148) 
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identity 
code 

one-cliaracter 
parity code" 

m-character 
Hamming code 

Hamming code 
and parity check* 
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upon {wciuu) 

wc = wd + 1 
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number of columns in the 
constructed sampling matrix'^ 

/ log2("ver) 

r7,Vlog2(dim Vo) 
f loga (river) 



log2('^ 



WC 

wc 
wc 

Wc 
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■ ujd + 'fn 

■ wu + m + l 



p = n 



log2(dim Vo) 


log2(n,) 


log2(n„e,) 


log2 (liver) 


Iog2(dim Vo) 


log2(n) 


log2(ll„e,) 


log2 (river) 



log2 



I = n V '°S2 (dim Vo ) log2 (n) 



^ log2(ridim Vo) 
^ Iog2(dimVn) 

/ 21og2(n) 
V log2(dim Vi 



)) 



y) 



Summary of symbol definitions 

p = number of sampling matrix columns 
n = number of sampling matrix rows 
nver = number of polytope vertices 
dim Vo = dimensionality of the polytope space 



Wc = length of codewords, thus n = (dimVQ)™° 
wu = length of data words, thus p = (nver)™" 
m = (integer) index of Hamming code 



" Parity characters are calculated mod river rather than the (more common in the literature) mod 2. 

Hamming-and-parity codes are sometimes called SECDED codes (single-error correct, double-error detect). 
■^The given functions p(n, river, dim Vo) are exact, and are valid for both real and complex vector spaces. 



Table 1: Recipes for deterministically constructing sampling matrices by the methods of 
Section The primary design variables are taken to be the dimensionality of polytope 

space dimVo, the number of polytope vertices n^er, and the desired number of sampling 
matrix rows n. All expressions are exact, and the results apply to both real and complex 
vector spaces. The expressions are organized so as to make manifest that increased minimal 
Hamming distance is associated with decreased column-length p; this is the central design 
trade-off. 



4.6.6. Petal-counting in GK geometry via coding theory These sampling theory results 
have a direct quantitative relation to the efflorescent GK geometry that we discussed in 
Sections 1.5.8 and 2.11 Specifically, we are now able to construct a petal-vector description 
of GK manifolds, and verify that they indeed have exponentially many petals. 

We consider a spin-^ rank-1 GK manifold having ngpin spins. The preceding tetrahedral 
GK construction deterministically generates a dictionary of petal- words {w^ : A; E 1 , 4"=?'" } 
in one-to-one correspondence with petal-vector states {\w^) : k G l,4"'=p'"}. This dictio- 
nary of states of course exponentially over-complete, since its number of words is 4"=p'" = 
dim'H. Yet we also know that random pairs of petal- vectors in our tetrahedral dictio- 
nary are pairwise orthogonal to an excellent approximation, because their median Hamming 
distance is 3nspin/4, and consequently from (152) their median pairwise wedge product is 

\w' ^w^\^ = 1 - 3-3"==p-/l 

As a concrete exercise in petal-counting, we consider a system of ngpin = 16 spin-^ 
particles. The tetrahedral construction generates a dictionary of 4^^ = 2^"^ petal-vectors for 
this system, each word of which labels a petal whose state- vector has a wedge separation of 
\w^' A w^\'^ > 2/3 from the state-vector of all other petals. A subset of that petal dictionary 
having minimal Hamming distance 4 is specified by the SECDED code of Table [TJ This 
SECDED subset has Hamming parameter m = 4, and hence m + 1 = 5 characters out of 
16 in each word are devoted to error-correcting. The resulting (smaller) error-corrected 
dictionary has 4^^~^ = 4^^ = 2^^ petal-vectors, and the sampling matrix whose columns are 
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the petal- vectors of this dictionary therefore has re = 2^^ rows and p 
column wedge products satisfy the (exact) pairwise inequality 



\w' Aw''\^ > 1 



for all i ^ k 



2^^ columns, whose 



(154) 



These calculations confirm our previous conclusion from Riemann curvature analysis, that 
even a rank-one GK manifold contains exponentially many petals. They also illustrate 
that the deterministic construction of high-quality sampling matrices involves sophisticated 
trade-offs in error-correcting codes. 

4-6.7. Constructing a Dantzig selector for quantum states We now have all the ingredients 
we need to establish the following remarkable principle: any quantum state that can be 
written as a (sparse) sum of petal-vectors can be recovered from sparse random projec- 
tions by convex programming methods. The method is as follows. We approximate the 
minimization problem ( 139[ ) that began this section in the form 



mm 



11^0 -^k(c) 



~ mm 

w 



\tpQ — Xw 



|2 
1/2 



(155) 



Here the columns of X are the petal-states of our dictionary (the preceding tetrahedral 
dictionary will do) and w is a column vector of petal coefficients (one coefficient for every 
vector in our dictionary). 

We then further approximate the above minimization problem in any of several standard 
forms [79], see also [35 ] [Ml [137] . These forms include (e.g.): 



min hi , s.t. 

w " "'1 



min i/jq 

w 



Xw\\, 

11(2 



s.t. 



x\-4^n- Xw)\\, <x 

< A 



w 



m.in\\w\\, +X\\iI)q — Xw\ 



Dantzig selector (156a) 
LASSO (156b) 
basis pursuit (156c) 



Here A is a parameter that is adjusted on a per-problem basis. Although the relative merits 
of the above optimizations are the subject of lively debate, for many practical problems 
they all work well. The Dantzig selector optimization (156a) in particular can be posed as 



an explicitly convex optimization problem that can be solved by a Dantzig-type simplex 
algorithm |37j (among other methods). 

To recapitulate, the key physical idea behind the first step (155) of the above two- 



step transformation is to represent a general state as a sparse superposition of petal-states. 



The key mathematical idea behind the second step (156a-c) is to approximate the result- 
ing sparse minimization problem as any of several forms that can be efficiently solved by 
numerical means. 

A second key mathematical idea is that the column dimensions of X can be very large — 
much larger than the Hilbert space dimension dimTC — provided that efficient algorithms 
exist for calculating the product Xw without calculating either X or iv explicitly (as was 
discussed earlier in Section 2.1). This is why deterministic methods for constructing X 



are essential to the feasibility of quantum optimization by Dantzig selection and related 
methods. 

This construction provides a non-trivial mathematical explanation of why the numerical 
optimizations of this article are well-behaved: the early coarse-grained, non-linear stages can 
be regarded as implicitly solving a convex optimization problem over petal-states, and the 
later fine-grained stages are solving a problem which is linear to a reasonable approximation. 

Boyd and Vandenberghe's textbook ambitiously asserts [211 in the Preface] 
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subsets with 


subsets with 


median A^in 


median Amin 


median A^ax 


median A^ax 


sparsity 


(s) 


Ss e (0, 1) 


Ss e (0, 1) 


of subsets 


of su bsets 


of subsets 


of subsets 


S 




(tetrahed ral)" 


(random) 


(tetrahedral) 


(random) 


(tetrahed ral) 


(random) 


1 


16 


16 (100%) 


16 (100%) 


1.000000 


1.000000 


1.000000 


1.000000 


2 


120 


120 (100%) 


120 (100%) 


0.666666 


0.691952 


1.333333 


1.308047 


3 


560 


560 (100%) 


545 (97%) 


0.490824 


0.445988 


1.509175 


1.594857 


4 


1820 


1804 (99%) 


1350 (74%) 


0.394501 


0.283165 


1.664513 


1.855300 


5 


4368 


3852 (88%) 


1376 (32%) 


0.265070 


0.173275 


1.815779 


2.104026 


6 


8008 


4408 (55%) 


413 (5%) 


0.168143 


0.096862 


1.971790 


2.341109 


7 


11440 


1360 (12%) 


21 (0%) 


0.100638 


0.044849 


2.103133 


2.572424 


8 


12870 


12 (0%) 


(0%) 


0.044634 


0.012520 


2.240730 


2.787555 


" Tetrahedral result 


s are exact (8 x 


16 sampling matrix, tetrahedral construction with Hamming distance 2). 


'Gaussian results are median values from 100 trials (8 x 16 random 


complex matrix with columns 


normalized). 



Table 2: RIP properties of 8 x 16 sampling matrices created via a deterministic tetrahedral 
construction, contrasted with same-size random Gaussian sampling matrices. By definition, 
matrices for which 100% of subsects have 6s G (0, 1) are RIP in order S. For all sparsities the 
tetrahedral construction yields RIP properties that are superior to random constructions. 



With only a bit of exaggeration, we can say that, if you formulate a practical problem 
as a convex optimization problem, then you have solved the original problem. 

But this assertion must be regarded with caution when it comes to convex optimization 
over quantum state-spaces, because the matrices and vectors involved are of enormously 
larger dimension than is usually the case in convex optimization. Figueiredo, Nowak, and 
Wright [To] and also Cai and Lv |29] discuss this domain, and it is clear that Cai and 
Lv's conclusion "Clearly, there is much work ahead of us" applies especially to compressive 
quantum sensing, sampling, and simulation. 

4-6.8. RIP properties of deterministic versus random sampling matrices It is clear from 
the preceding discussion that over-complete dictionaries of word-states having the 

approximate orthogonality property {w^\w^) = {XX^)ij ~ 6ij are desirable both for simula- 
tion purposes and for sampling purposes. Stimulated by the work of Candes and Tao [37] . 
an extensive and rapidly growing body of work characterizes such matrices in terms of the 
restricted isometry property (RIP). We now briefly discuss the RIP of tetrahedral sampling 
matrices, mainly following the notation and discussion of Baraniuk et al. [lOj . 

We regard the word indices i and j in {w'^\w^) as the rows index and column index of 
a Hermitian matrix. We specify a subset T of word indices, and we define the sparsity S 
of that subset to be 5 = #r. Then {w'^\w^)t = {w^lw-') : i,j £ T is an S x S Hermitian 
matrix, which we take to have minimal (maximal) eigenvalues Amin (Amax)- Then the 
isometry constant 6s of Candes and Tao is by definition 

6s = max(l - Amin, Amax " 1) (157) 

Our word-state dictionary is said to have the restricted isometry property for order S iff 
6s S (0, 1) for all subsets T having sparsity S. Physically speaking, a dictionary of p word- 
states having the RIP property in order S has the property that any set of S words is 
(approximately) mutually orthogonal. 

Testing for the RIP property is computationally inefficient, since (at present) no known 
algorithm is significantly faster than directly evaluating Amin and Amax for all (^) distinct 
subsets T. Referring to Table [ij we see that a spin-^ tetrahedral dictionary of three-letter 
words, one of which is a parity-check character, such that the minimal Hamming distance 
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is two, yields a sampling matrix having p = 4p = 16 columns and with n = 2^ = 8 
rows. To calculated the RIP properties of this dictionary, the maximum sparsity we need to 
investigate is 5 = n = 8, for which (^) = ( g^) = 12870 subsets must be evaluated, which is a 
feasible number. As summarized in Table |2] the tetrahedral construction yielded sampling 
matrices having the RIP property for sparsity S = 1, 2, 3, while for higher values of S the 
fraction of subsets having 6s G (0, 1) dropped sharply. 

For purposes of comparison, we computed also the median RIP properties of 8 x 16 ran- 
dom Gaussian matrices. We found for all values of sparsity, the RIP properties of Gaussian 
random matrices were inferior to those of the deterministic tetrahedral construction. We are 
not aware of any previous such random-versus-deterministic comparisons in the literature. 
Since is is known that the Gaussian random matrices are RIP in the large=p limit, we were 
surprised to find that their RIP properties are unimpressive for moderate values of p. 

In preliminary studies of larger matrices, we found that known asymptotic expressions 
for the extremal singular values of Gaussian random sampling matrices — due to Marcenko 
and Pastur |134j . Geman [57], and Silverstein |183j . as summarized for CS purposes by 
Candes and Tao |33| see their Sec. Ill] — were empirically accurate for tetrahedral sampling 
matrices too, for all values of the row dimension n < 256 and all values of the sampling 
parameter S < n. 

We emphasize however that although the average-case performance of these petal-vector 
sampling matrices is empirically comparable to Gaussian sampling matrices, their worst- 
case performance is presently unknown, and in particular such key parameters as their 
worst-case isometry constants are not known. 

As Baraniuk et al. [lOj note: "the question now before us is how can we construct matri- 
ces that satisfy the RIP for the largest possible range of 5." It is clear that answering this 
question, in the context of the deterministic geometric construction given here, comprises 
a challenging problem in coding theory, packing theory, and spectral theory, involving so- 
phisticated trade-offs among the competing goals of determinate construction, large (and 
adjustable) p/n ratio in the sampling matrix, and small isometry constants for all values of 
the sparsity parameter S <n. 

4-6.9. Why do CS principles work in QMOR simulations? Guided by the preceding analy- 
sis, we now try to appreciate more broadly why CS principles "work' in QMOR simulations 
by systematically noting mathematical parallels between the two disciplines. We will see 
that these parallels amount to an outline for extending the mathematical foundations of CS 
to provide foundations for QMOR-CS. 

As our first parallel, we remark that what Candes and Tao call ^37j compressible ob- 
jects are ubiquitous in both the classical and quantum worlds. This ubiquity is not easily 
explained classically, and so almost always it is simply accepted as a fact of nature; for 
example almost all visual fields of interest to human beings are compressible images. In 
contrast, we have seen that the ubiquity of quantum compressible objects has a reasonably 
simple explanation: most real-world quantum systems are noisy, and noisy systems can 
be modeled as synoptic measurement processes that compress state trajectories; working 
through the mathematical details of this synoptic compression was of course our main con- 
cern in Section |3j From this quantum informatic point of view, it is a fundamental law 
of nature that any quantum system that has been in contact with a thermal reservoir (or 
equivalently, a measurement-and-control system) is a compressible object. 

The second parallel is the availability of what the CS field calls dictionaries |19Q) of 
the natural elements onto which both classical and quantum compressible objects are pro- 
jected. For example, wavelet dictionaries are well-suited to image reconstruction. In the 
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QMOR formalism of this article, the parallel quantum dictionary is (of course) the class of 
multilinear biholomorphic GK polynomials that define the Kahlerian geometry of QMOR 
state-spaces (Section [2]). This is not a linear dictionary of the type generally discussed 
in the CS literature, but rather is an algebraic generalization of such dictionaries. In the 
language of Donoho [61], open quantum systems exhibit a generalized transform sparsity 
whose working definition is the existence of high-fidelity projections onto GK manifolds. 

The third parallel is the existence of robust, numerically efficient methods for projec- 
tion and reconstruction. It is here that the mathematical challenges of aligning QMOR 
with CS are greatest. In our own research we have tried non-CS/non-QMOR optimization 
techniques — like regarding ■j/>q — ip(c) = as the definition of an algebraic variety, and 
decomposing it into a Groebner basis — but in our hands these methods perform poorly. 
Turning this observation around, it is possible that the efficient methods of QMOR-CS 
might find application in the calculation of (specialized algebraic forms of) Groebner bases. 

Although a substantial body of literature exists [189] for minimizing functions that are 
convex along geodesic paths on Riemannian state-spaces — which generalizes the notion of 
convexity on Euclidean spaces — there does not seem to be any similar body of literature on 
the convexity properties of holomorphic functions on Kahlerian state-spaces. 

We have previously quoted Shing Tung Yau's remark |200| p. 21]: "While we see great 
accomplishments for Kahler manifolds with positive curvature, very little is known for 
Kahler manifolds [having] strongly negative curvature." By the preceding construction, 
we now appreciate that (negatively curved) GK manifolds have embedded within them 
lattices that display all the intricate mathematical structure of coding theory — so that it 
is not surprising that the geometric properties of these manifolds resists easy analysis. It 
seems that the ultrahigh-dimensional model selection of Cai and Lv [29] can be described — 
with more-or-less equal mathematical justification — in terms of the differential geometry 
of ruled manifolds, or alternatively in terms of coding theory, or alternatively in terms of 
optimization theory. 

There is also the as-yet unexplored practical issue of whether quantum optimization ofl2- 



type functions over GK polynomials like (139) is more efficient, less efficient, or comparably 



efficient to CS-type optimization over petal- words of /i-type functions like (156a-c). This 
question is analogous to the long-standing issue in CS of whether interior-point methods 
are superior to edge- and- vertex polytope methods ... the answer after five decades of CS 
research being "yes, sometimes." 



Efficient numerical means for evaluating the Penrose pseudo-inverse of ( 143 1 are needed, 
as this inversion is the most computationally costly step of our sparse reconstruction codes 
as they are presently implemented. Preconditioned conjugate gradient techniques are one 
attractive possibility |5H IMl [96] , because these techniques lend themselves well to the large- 
scale parallel processing. The algebraic structure of the GK metric tensor creates additional 
algorithmic challenges and opportunities that (so far as the authors are aware) have not 
been addressed in the computing literature. 

Finally, suites of test problems and open-source software tools have contributed greatly 
to the rapid development of CS theory and practice |190] . and it would be valuable to have 
a similar suite of problems and tools for the simulation of open quantum systems. 

5. Conclusions 

As Terence Tao has remarked |187] 

The field of [partial differential equations] has proven to be the type of mathematics 
where progress generally starts in the concrete and then flows to the abstract, rather 
than vice versa. 



5.1 Concrete applications of large-scale quantum simulation 



91 



The recipes in this article have resulted from a flow in the opposite direction, from abstract 
to concrete, in which abstract ideas from quantum information theory, algebraic geometry, 
quantum physics, and compressive sensing have found find concrete embodiment in practical 
recipes for large-scale quantum simulation. 

Key abstract ideas from quantum information theory (QIT) Our recipes have adopted 
from quantum information theory the key idea that noise processes can be modeled as 
covert measurement processes. This leads naturally to the idea that quantum states that 
have been in contact with a thermal reservoir (or equivalently, a measurement and control 
process) are compressible objects. 

Key ideas from algebraic geometry Our recipes have adopted from algebraic geometry the 
key idea that reduced-order quantum state-spaces can be described as geometric objects, 
using the language and methods of algebraic and differential geometry. In particular, quan- 
tum trajectories can be described in terms of drift and diffusion processes upon state-space 
manifolds, just as in classical modeling and simulation theory. 

Key ideas from quantum physics theory Our recipes have adopted from theoretical quan- 
tum physics the key idea that quantum states have multiple unravelings, and that the 
efficiency of a calculation can be optimized by choosing an appropriate unravelling. 

Key ideas from quantum physics experiments Our recipes have adopted from experimental 
quantum physics the key idea that mathematical ingredients of quantum simulation map 
one-to-one onto familiar physical systems such as measuring devices, and also the opera- 
tional principle that the main deliverable of a quantum simulation is accurate prediction of 
the results of physical measurements. 

Key ideas from compressive sensing, sampling, and simulation ( CS) Our recipes have 
adopted from CS the idea that optimization problems involving compressible objects (like 
quantum states) can often be transformed into convex optimization problems. This can 
lead to both faster algorithms and improved physical insight. 

5.1. Concrete applications of large-scale quantum simulation 

By combining the preceding abstract ideas, the objective that began this article 

... to enable the reader to design and implement practical quantum simulations, guided 
by an appreciation of the geometric, informatic, and algebraic principles that govern 
simulation accuracy, robustness, and efficiency. 

has now been achieved in a preliminary sense, albeit there is much further work to be done. 
Now we consider some practical applications. 

5.1.1. The goal of atomic-resolution biomicroscopy The goal of atomic-resolution micro- 
scopy was the main motivation for developing the simulation algorithms described in this 
article. This goal was proposed as early as 1946 by Linus Pauling who envisioned "If it were 
possible to make visible the individual molecules of the serum proteins and other proteins of 
similar molecular weight, all the uncertainty which now exists regarding the shapes of these 
molecules would be dispelled" |151] . Later that same year, John von Neumann (possibly 
speaking as a reviewer of Pauling's proposal |116] ) wrote a letter to Norbert Wiener [145] 
that embraced and extended Pauling's vision. The letter expressed a strikingly modern 
vision of atomic-level structural and systems biology: 

There is no telling what really advanced electron-microscopic techniques will do. . . . 
A "true" understanding of [viral-scale] organisms may be the first step forward and 
possibly the greatest step that may at all be required. I would, however, put on "true" 
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understanding the most stringent interpretation possible: That is, understanding the 
organism in the exacting sense in which one may want to understand a detailed drawing 
of a machine, i.e. finding out where every individual nut and bolt is located. 

It was not until 1959 that Richard Feynman — who spent a sabbatical year working as 
a biochemist [63] — issued his famous challenge: "Is there no way to make the electron 
microscope more powerful? . . . Make the microscope one hundred times more powerful, 
and many problems of biology would be made very much easier" ^68j . 

Unfortunately, the problem of electron-beam radiation damage to fragile biological mol- 
ecules proved intractable |11U| . and so the Pauling- von Neumann- Feynman challenge of 
achieving atomic-resolution biomicroscopy remained unanswered for several decades. 

New ideas were needed, and three key ideas that emerged in ensuing decades were mag- 
netic resonance imaging, nanotechnology, and quantum measurement theory. The early 
stages of development of each of these new fields was slow, not because fundamentally new 
concepts of mathematics or physics were required — the key concepts were reasonably famil- 
iar to Pauling, von Neumann, and Feynman's generation — but because each field sought 
to push familiar concepts to extreme limits. Magnetic resonance was a familiar concept; 
exploiting the tiny magnetic resonance signals for 3D imaging purposes was novel. Making 
devices smaller was a familiar concept; fabricating micron-scale and nanometer-scale devices 
was novel. Quantum measurement was a familiar concept; studying in detail the quantum 
evolution of small, continuously observed systems was novel. 

These research fields are united in magnetic resonance force microscopy (MRFM) , which 
was conceived explicitly as a means of meeting the Pauling-von Neumann-Feynman chal- 
lenge of atomic-resolution biomicroscopy jl77| I178| I179| [181^ . 

5.2. The acceleration of classical and quantum simulation capability 

Simulation technologies, both classical and quantum, began an immense surge of progress 
during the Pauling-von Neumann-Feynman era, and this surge has continued to the present 
day. This is true especially in the classical domain, where simulation tools have become 
essential to system- level engineering \20\ 1113] . In an era in which a new aircraft, a new 
processor chip, or a new drug can readily incur system development costs in excess of one 
billion dollars, engineering processes that maximize confidence, reliability, and economy 
have become a practical necessity. 

Feynman, in his seminal 1982 article Simulating physics with computers [67], argued 
that generic quantum physics problems are exponentially hard to simulate on a classical 
computer. But Feynman's analysis was vague about what constitutes generic quantum 
physics. The viewpoint we have developed in this article is that any quantum system that 
has been in contact with a thermal reservoir is — in principle at least — a compressible object 
(in the language of CS theory) and thus is amenable to simulation with classical resources. 

The present status of quantum simulation algorithms has, in our experience, striking 
parallels to the status of linear programming algorithms during 1947-2004 [561. It was clear 
for many decades that linear programming methods were exceedingly useful for solving prac- 
tical problems, but their mathematical foundations in convex set theory were established 
only very slowly (see Spielman and Teng's 2004 article |185j for a review and a reasonably 
definitive solution). 

This slow-but-steady progress illustrates Dantzig's principle [S]: 

In brief, one's intuition in higher dimensional space is not worth a damn! Only now, 
almost forty years after the time the simplex method was first proposed, are people 
beginning to get some insight into why it works as well as it does. 
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and also Feynman's words [66] 

I think the problem is not to find the best or most efficient method to proceed to a 
discovery, but to find any method at all. . . . [That is why] it is useful to have a wide 
range of physical viewpoints and mathematical expressions of the same theory. 

5.3. The practical realities of quantum system engineering in MRFM 

The practical experience of operating an MRFM devices is very much like operating a small 
satellite that is distant from the experimenter not in space, but in scale. In particular, the 
state-space of the MRFM device includes the quantum state-space of spins in the sample. 
The QMOR analysis methods of this article were conceived specifically to allow the efficient 
a priori modeling of this quantum state-space, by methods that have the potential to sub- 
stantially extend present capabilities in modeling large-scale spin-systems To the extent 
that future progress in QMOR analysis allows this goal to be achieved, then the optimiza- 
tion of MRFM technology can proceed partially in silico, which will help retire technical 
risk, speed development, and build alliances among sponsors, researchers, enterprises, and 
customers. 

5.4. Future roles for large-scale quantum simulation 

In our view, the single most important role for quantum system engineering (QSE) will be to 
help sustain the exponentially cumulative technological progress that characterized the 20st 
century. In the context of computing this exponentiation is popularly known as "Moore's 
Law," but it is fair to say that similar exponentially cumulative progress is evident in fields 
such as nanotechnology, information theory, and (especially) biology. 

A recent theme issue of IBM Journal of Research and Development describes large- 
scale simulation codes running on "Blue Gene" hardware that is approaching petafiop-scale 
computation speeds [HI- Both classical [Ml [721 [731 (121 [SS] and quantum [HI ^ [TUH] 
simulations are reviewed, and it is fair to say that the boundary between these two kinds of 
simulations is becoming indistinct, in particular when it comes to computing inter- atomic 
potentials that are both numerically efficient (classical) and accurate (quantum). This 
continues a sixty-year record of mutually supportive progress in hardware, software, and 
algorithm development |60j . 

From a geometric point of view, modern multi-processor computer architectures are ex- 
ceptionally well-suited to the efficient computation of fundamental geometric objects such as 
polytopes, metrics, drift vectors, gradients, and diffusion tensors. These fundamental geo- 
metric objects are the raw building blocks — both conceptually and as software libraries — for 
broad classes of system simulations. In particular, the recipes of this article demonstrate 
that both classical and quantum systems can be simulated using a shared set of abstract 
mathematical ideas and concrete software tools that are well-matched to distributed com- 
puting architectures. 

In summary, twenty-first century technologies seek to maximize the pace, coordination, 
and reliability of technology development, to create products that press against the quantum 
quantum and thermodynamic limits of device speed, sensitivity, size and power. For all such 
technologies, the quantum simulation recipes of this article promise to be useful. 

Acknowledgements 

Author J. A. Sidles gratefully acknowledges the insights gained from discussions with Al 
Matsen on quantum chemistry. Rick Matsen on regenerative medicine, and Eric Matsen on 
informatic biology. 

This work was supported by the Army Research Office (ARO) MURI program under 



5.4 Future roles for large-scale quantum simulation 



94 



contract #W911NF-05-l-0403, and by the National Institutes of Health NIH/NIBIB under 
grant #EB02018. Prior support was provided by the ARO under contact #DAAD19-02-l- 
0344, by IBM under subcontract #A0550287 through ARO contract #W911NF-04-C-0052 
from DARPA, and by the NSF/ENG/ECS under grant #0097544. 



REFERENCES 



95 



References 

[1] Abramowitz M and Stegun I A, eds 1964 Handbook of Mathematical Functions, 10th 
Edition, US Government Printing Office Washington. 

[2] Adler S L 2000 Physics Letters A pp. 58-61. 

[3] Aharonov D and Ben-Or M 1996 37Th Annual IEEE Symposium on Foundations of 
Computer Science (FOCS'96) pp. 46-55. 

[4] Al-Hassanieh K, Dobrovitski V, Dagotto E and Harmon B 2006/07/21 Phys. Rev. Lett. 
(USA) 97(3), 037204- 1. 

[5] Albers D J and Reid C 1986 The College Mathematics Journal 17(4), 292-314. 

[6] AHcki R and Lendi K 1987 Quantum Dynamical Semigroups and Applications Springer- 
Verlag. 

[7] Antoulas A C 2005 Approximation of Large-Scale Dynamical Systems SIAM. 

[8] Apostolov V and Gualtieri M 2007 Communications in Mathematical Physics 
271(2), 561-575. 

[9] Ashtekar A and SchiUing T A 1999 in 'On Einstein's Path' Springer pp. 23-65. 
arXiv : gr-qc/9706069. 

[10] Baraniuk R, Davenport M, DeVore R and Wakin M 2008 A simple proof of the re- 
stricted isometry property for random matrices, to appear in Constructive Approxi- 
mation. 

[11] Baraniuk R G and Wakin M B 2008 Random projections of smooth manifolds, to 
appear in Foundations of Computational Mathematics. 

[12] Belcastro S M and Yackel C 2007 Making Mathematics with Needlework A. K. Peters 
Ltd. 

[13] Benvegnu A, Sansonetto N and Spera M 2004 Journal of Geometry and Physics 
(Netherlands) 51(2), 229-43. 

[14] Besse A L 1987 Einstein Manifolds Springer- Verlag. 

[15] Beylkin G and Mohlenkamp M J 2002 Proceedings of the National Academy of Science 
99(16), 10246-10251. 

[16] Beylkin G and Mohlenkamp M J 2005 SIAM Journal on Scientific Computing 
26(6), 2133-2159. 

[17] Bickcl P J 2007 Ann. Statist. 35(6), 2352-2357. 

[18] Blatt R and Zoller P 1988 European Journal of Physics 9(4), 250-256. 

[19] Bohm E, Bhatele A, Kale L V, Tuckerman M E, Kumar S, Gunnels J A and Martyna 
G J 2008 IBM Journal of Research and Development 52(1/2), 159-176. 

[20] Booton R and Ramo S 1984 IEEE Transactions on Aerospace and Electronic Systems 
AES-20, 306-9. 



REFERENCES 96 



[21 

[22 

[23; 
[24; 
[25 

[26; 

[27; 

[28; 
[29; 
[30 
[31 

[32 
[33 
[34; 
[35 

[36 
[37; 

[38; 
[39; 

[4o; 

[41 
[42 



Boyd S and Vandenberghe L 2004 Convex Optimization Cambridge University Press 
Cambridge. 

Braginsky V R and Khalili F Y 1992 Quantum Measurement Cambridge University 
Press. 

Braginsky V, Vorontsov Y and Thorne K 1980 Science (USA) 209(4456), 547 - 57. 
Bravyi S and Kitaev A 2005 Physical Review A 71, 022316. 

Breuer H P and Petruccione F 2002 The Theory of Open Quantum Systems Oxford 
University Press. 

Bruland K J, Garbini J L, Dougherty W M, Chao S H, Jensen S E and Sidles J A 1999 
Review of Scientific Instruments 70, 3542-3544. 

Buhrman H, Cleve R, Laurent M, Linden N, Schrijver A and Unger F 2006 ^7T/i Annual 
IEEE Symposium on Foundations of Computer Science (FOCS'06) pp. 411-419. 

Burke K, Werschnik J and Gross E K U 2005 Journal of Chemical Physics 123, 062206. 

Cai T T and Lv J 2007 Ann. Statist. 35(6), 2365-2369. 

Candes E J and Romberg J 2006 Found. Comput. Math. 6(2), 227-254. 

Candes E J, Romberg J K and Tao T 2006a Comm. Pure Appl. Math. 59(8), 1207- 
1223. 

Candes E J, Romberg J and Tao T 20066 IEEE Trans. Inform. Theory 52(2), 489-509. 

Candes E J and Tao T 2005 IEEE Trans. Inform. Theory 51(12), 4203-4215. 

Candes E J and Tao T 2006 IEEE Trans. Inform. Theory 52(12), 5406-5425. 

Candes E and Justin Romberg C 2005 Zi-magic : Recovery of sparse signals Technical 
report California Institute of Technology. Software manual, available at http : / / www .] 
|acm. caltech. edu/llmagic/j 

Candes E and Romberg J 2007 Inverse Problems 23(3), 969-985. 
Candes E and Tao T 2007a Annals of Statistics 35(6), 2313-2351. 
Candes E and Tao T 20076 Ann. Statist. 35(6), 2392-2404. 



Capelle K 2006 Brazilian Journal of Physics 36(4A), 1318-43. arXiv: cond-mat/ 
I02ll443v5l 

Carmichael H J 1993 An Open Systems Approach to Quantum Optics Lecture Notes 
in Physics, New Series M18 Springer- Ver lag. 

Carmichael H J 1999 Statistical Methods in Quantum Optics I: Master Equations and 
Fokker-Planck Equations Springer. 

Cavalcante R G 2002 Performance Analysis of Signal Constellations in Riemannian 
Varieties PhD thesis FEEC-UNICAMP, Brazil. 



REFERENCES 



97 



[43] Cavalcante R G and Jr. R P 2005 Computational & Applied Mathematics 24(2), 173- 
92. 

[44] Cavalcante R G, Lazari H, Lima J D and Jr. R P 2005 in 'Algebraic Coding Theory 
and Information Theory (DIMACS Series in Discrete Mathematics and Theoretical 
Computer Science Vol. 68)' pp. 145-77. 

[45] Caves C 1982/10/15 Phys. Rev. D, Part. Fields (USA) 26(8), 1817 - 39. 

[46] Caves C M 2000 Completely positive maps, positive maps, and the Lindblad form. 
On-line memorandum to C. A. Fuchs and J. Renes: ^http : / / inf o . phys . unm . 
]2^caves/ report s/lindblad . pdfj 

[47] Chinmay Hegde M B W and Baraniuk R G 2007 Random projections for manifold 
learning: Proofs and analysis Technical Report TREE 0710 Rice University, Depart- 
ment of Electrical and Computer Engineering. 

[48] Choi M D 1975 Linear Algebra and Rs Applications 10, 285-290. 

[49] Chow B 2004 The Ricci Flow : An Introduction American Mathematical Society. 

[50] Chow B 2007 The Ricci Flow : Techniques and Applications American Mathematical 
Society. 

[51] Ciegis R 2005// Informatica (Lithuania) 16(3), 317 - 32. 

[52] Corless R M, Gonnet G H, Hare DEC, Jeffrey D J and Knuth D E 1996 Adv. Comput. 
Math. 5(4), 329-359. 

[53] Cox D, Little J and O'Shea D 2007 Ideals, Varieties, and Algorithms Undergraduate 
Texts in Mathematics, Third Edition, Springer New York. 

[54] Cramer C J 2004 Essentials of Computational Chemistry, 2nd Edition, Wiley. 

[55] Daley A J, KoUath C, Schollwock U and Vidal G 2004 Journal of Statistical Mechanics: 
Theory and Experiment pp. P04005-I-28. 

[56] Dantzig G B 1963 Linear programming and extensions Princeton University Press. 

[57] Davies E B 1976 Quantum Theory of Open Systems Academic Press. 

[58] del Castillo GET 2003 3-D Spinors, Spin-weighted Functions, and their Applications 
Birkhauser. 

[59] Dirac P A M 1930 Proceedings of the Cambridge Philosophical Society 26, 376-385. 

[60] Dongarra J J 2004 'An interview with jack j. dongarra' Interview conducted by the 
Society for Industrial and Applied Mathematics, by Thomas Haigh. 

[61] Donoho D L 2006 IEEE Transactions on Information Theory 52(4), 1289 - 1306. 

[62] Donoho D and Stodden V 2006 in 'Neural Networks, 2006' International Joint Confer- 
ence on Neural Networks pp. 1916-1921. 

[63] Edgar R S, Feynman R P, Klein S, Lielausis I and Steinberg C M 1962 Cenetics 
47(2), 179-186. 



REFERENCES 98 

[64] Efron B, Hastie T and Tibshirani R 2007 Ann. Statist. 35(6), 2358-2364. 

[65] Ethier S, Tang W M, Walkup R and Oliker L 2008 IBM Journal of Research and 
Development 52{l/2), 105-116. 

[66] Feynman R P 1966 Physics Today p. 31. 

[67] Feynman R P 1982 International Journal of Theoretical Physics 21(6-7), 467-488. 

[68] Feynman R P 1992 Journal of Microelectromechanical Systems 1(1), 60-66. Transcript 
of lectm'e given at the December 1959 American Physical Society annual meeting. 

[69] Feynman R P, Leighton R B and Sands M 1965 The Feynman Lectures on Physics III: 
Quantum Mechanics Addison- Wesley. 

[70] Figueiredo M, Nowak R and Wright S 2007/12/ IEEE J. Sel. Top. Signal Process. 
(USA) 1(4), 586 - 97. 

[71] Fischer G and Piontkowski J 2001 Ruled Varieties: An Introduction to Algebraic Dif- 
ferential Geometry Vieweg & Sohn. 

[72] Fisher R T, Kadanoff L P, Lamb D Q, Dubey A, Plewa T, Calder A, Cattaneo F, 
Constantin P, Foster I, Papka M E, Abarzhi S I, Asida S M, Rich P M, Glendenin 
C C, Antypas K, Sheeler D J, Reid L B, Gahagher B and Needham S G 2008 IBM 
Journal of Research and Development 52(1/2), 127-136. 

[73] Fitch B G, Rayshubskiy A, Eleftheriou M, Ward T J C, Giampapa M E, Pitman 
M C, Pitera J W, Swope W C and Germain R S 2008 IBM Journal of Research and 
Development 52{l/2), 145-158. 

[74] Flaherty E J 1976 Hermitian and Kdhlerian Geometry in Relativity Springer- Ver lag. 

[75] Flanders H 1963 Differential Forms with Applications in the Physical Sciences Aca- 
demic Press. 

[76] Frenkel J 1934 Wave Mechanics, Advanced General Theory Clarendon Press, Oxford. 

[77] Frenkel T 2004 The Geometry of Physics: an Introduction (Second Edition) Cambridge 
University Press. 

[78] Friedan D H 1985 Annals of Physics 163, 318-419. 

[79] Friedlander M P and Saunders M A 2007 Annals of Statistic 35(6), 2385-2391. 

[80] Friesner R A 2005 Proceedings of the National Academy of Sciences 102(19), 6648-53. 

[81] G. Obinata B D O A 2000 Model Reduction for Gontrol System Design Springer. 

[82] Gallot S, Hulin D and Lafontaine J 2004 Riemannian Geometry Springer. 

[83] Garcia-Mata I, Saraceno M, Spina M E and Carlo G 2000 Phase space contraction and 
quantum operations, eprint: lquant-ph/0508190l 

[84] Gardiner C W 1985 Handbook of Stochastic Processes 2nd Edition, Springer. 

[85] Gardiner C W and ZoUer P 2000 Quantum Noise 2nd Edition, Springer. 



REFERENCES 99 

[86] Gauss C F 1827 Commentationes Societatis Regiae Scientiarum Gottingensis pp. 99- 
146. 

[87] Geman S 1980 Ann. Probab. 8(2), 252-261. 

[88] Gerlach W and Stern O 1921 Zeitschrift fur Physik 8(2), 110 - 111. 

[89] Gerlach W and Stern O 1922a Zeitschrift fur Physik 9(6), 349 - 352. 

[90] Gerlach W and Stern O 19226 Zeitschrift fur Physik 9(6), 353 - 355. 

[91] Goldberg S I and Kobayashi S 1967 Journal of Differential Geometry 1, 225-233. 

[92] Gottesman D 1997 Stabilizer Codes and Quantum Error Correction PhD thesis Cali- 
fornia Institute of Technology. 

[93] Gottfried K 1966 Quantum Mechanics W. A. Benjamin. 

[94] Gravvanis G and Giannoutakis K 2005// Sci. Program. (Netherlands) 13(2), 79 - 91. 

[95] Griffiths D J 2005 Introduction to Quantum Mechanics Pearson Prentice Hall. 

[96] Gupta A, Kumar V and Sameh A 1995/05/ IEEE Trans. Parallel Distrib. Syst. (USA) 
6(5), 455 - 69. 

[97] Gygi F 2008 IBM Journal of Research and Development 52(1/2), 137-144. 

[98] Hairer E, Lubich C and Wanner G 2006 Geometric Numerical Integration: Structure- 
Preserving Algorithms for Ordinary Differential Equations 2nd Edition, Springer. 

[99] Hardy L 1998 Contemporary Physics (UK) 39(6), 419 - 29. 

[100] Harris J 1992 Algebraic Geometry: A First Course Springer- Verlag. 

[101] Harrow A and Nielsen M 2003 Physical Review A - Atomic, Molecular, and Optical 
Physics 68, 012308. 

[102] Hartree D R 1928 Proceedings of the Cambridge Philosophical Society 24, 89-132. 

[103] Haus H A and Mullen J A 1983 in J. A Wheeler and W. H Zurek, eds, 'Quan- 
tum Theory and Measurement' Princeton series in physics Princeton University Press 
Princeton, New Jersey pp. 736-742. first published as |104] . 

[104] Haus H and Mullen J 1962 Physical Review 128(5), 2407- 2413. collected in [T05] . 

[105] Hecht J 2005 Laser Focus World (USA) 41(12), 71 - 4. 

[106] Heffner H 1962 Proceedings of the Institute of Radio Engineers (IRE) 50(7), 1604- 
1608. collected in |107j. 

[107] Heffner H 1983 in J. A Wheeler and W. H Zurek, eds, 'Quantum Theory and Mea- 
surement' Princeton series in physics Princeton University Press Princeton, New Jersey 
pp. 725-735. first published as [106j. 

[108] Helmer J and MuUer M 1958/04/ Institute of Radio Engineers Transactions on Mi- 
crowave Theory and Techniques MTT-6(2), 210 - 214. 



REFERENCES 



100 



[109] Henderson D and Taimina D 2001 Mathematical Intelligencer 23(2), 17-18. 

[110] Henderson R 1995 Quarterly Reviews of Biophysics 28(2), 171-193. 

[Ill] Hou B Y and Ho B Y 1997 Differential Geometry for Physicists World Scientific. 

[112] IBM Blue Gene Team T 2008 IBM Journal of Research and Development 
52(1/2), 199-220. 

[113] Johnson F T, Tinoco E N and Yu N J 2005/12/ Computers and Fluids 34(10), 1115- 
51. 

[114] Jones R O and Gunnarsson O 1989 Reviews of Modem Physics 61(3), 689-746. 

[115] Jost J 1997 Nonpositive Curvature: Ceometric and Analytic Aspects Biikhausei . 
[116] Kay L 1998. Personal communication. 

[117] Kay L E 1993 The Molecular Vision of Life: Caltech, the Rockefeller Foundation, and 
the Rise of the New Biology Oxford University Press. 

[118] Klauder C and Skagerstam B S 1984 in 'Coherent States' World Scientific Publishing 
pp. 3-74. 

[119] Klumper A, Schadschneider A and Zittartz J 1993/11/01 Europhys. Lett. (Switzer- 
land) 24(4), 293 - 7. 

[120] Kobayashi S 1987 Differential Geometry of Complex Vector Bundles Iwanami Shoten, 
and Princeton University Prcsss. 

[121] Kohn W 1999 Reviews of Modern Physics 71(5), 1253-66. 

[122] Kraus K 1971 Annals of Physics 64, 311-335. 

[123] Kraus K 1983 States, Effects, and Operations: Fundamental Notions of Quantum 
Theory Springer. 

[124] Kumar S, Huang C, Zheng G, Bohm E, Bhatele A, Phillips J C, Yu H, and Kale L V 
2008 IBM Journal of Research and Development 52(1/2), 177-188. 

[125] Leggett A 2002/04/22 J Phys., Condens. Matter (UK) 14(15), 415 - 51. 

[126] Leggett A 2002a Phys. Scr. Vol. T (Sweden) T102, 69 - 73. 

[127] Leggett A 20026 in 'Fundamentals of Quantum Information: Quantum Computation, 
Communication, Decoherence and All That' Johannesburg, South Africa pp. 3 - 45. 

[128] Leggett A J 2004 in A. J Leggett, B Ruggiero and P Silvestrini, eds, 'Quantum 
Computing and Quantum Bits in Mesoscopic Systems' Third International Workshop 
on Macroscopic Quantum Coherence and Computing (Naples 2002) Kluwer Academic 
pp. 1-12. 

[129] Leggett A J, Chakravarty S, Dorsey A T, Fisher P A, Garg A and Zwerger W 1987 
Rev. Mod. Phys. (USA) 59(1), 1 - 85. 

[130] Leibfried D, Blatt R, Monroe C and Wineland D J 2003 Reviews of Modem Physics 
75(1), 281-324. 



REFERENCES 101 

[131] Lindblad G 1976 Communications in Mathematical Physics 48, 119-130. 

[132] Loan C F V and Golub G H 1996 Matrix Computations, Third Edition, Johns Hopkins 
University Press. 

[133] Lubich C 2004 Mathematics of Computation 74(250), 765-779. 

[134] Marcenko V A and Pastur L A 1967 Matematicheskiui Shornik. Novaya Seriya 72 
(114), 507-536. 

[135] Martin D 2002 Manifold Theory: An Introduction for Mathematical Physicists Hor- 
wood Pubhshing Limited. 

[136] McLachlan A D 1964 Molecular Physics 8, 39-44. 

[137] Meinshausen N, Rocha G and Yu B 2007 Ann. Statist. 35(6), 2373-2384. 

[138] Mermin N 1985 Phys. Today (USA) 38(4), 38 - 47. 

[139] Messiah A 1966 Quantum Mechanics, Volume //John Wiley & Sons. 

[140] Mohlenkamp M J and Monzon L 2005 Mathematical Intelligencer 27(2), 65-69. 

[141] Moore S 2007 IEEE Spectrum (USA) 44{3), 15 - 17. 

[142] Moroianu A 2007 Lectures on Kdhler Ceometry Cambridge University Press. 
[143] Morris N F 1977 Computers & Structures 7(1), 65-72. 

[144] MuUe U and van de Ven A E M 1999 Ceneral Relativity and Cravitation 31(11), 1759- 
68. 

[145] Neumann J V 1997 in V Mandrekar and P. R Masani, eds, 'Proceedings of the Nor- 
bert Wiener Centenary Congress, 1994' Vol. 52 of Proceedings of Symposia in Applied 
Mathematics American Mathematical Society, Providence, RI pp. 506-512. 

[146] Nielsen M A and Chuang I L 2000 Quantum Computation and Quantum Information 
Cambridge Press. 

[147] Noor A K 1982 Computer Methods in Applied Mechanics and Engineering 34((1- 
3)), 955-985. 

[148] Noor A K 1994 Applied Mechanics Reviews 47(5), 125-145. 

[149] on Mathematical Challenges From Computational Chemistry C 1995 Mathematical 
Challenges from Theoretical Quantum Chemistry National Academy Press, National 
Research Council. 

[150] Papenbrock T, Juodagalvis A and Dean D J 2004 Physical Review C69(2), 024312. 

[151] Pauling L 1946 'The possibilities for progress in the fields of biology and biologi- 
cal chemistry'. Proposal to Warren Weaver, Rockefeller Foundation, June 19, 1946. 
This file was supplied by historian Lily Kay |117j : the original sources are avail- 
able through the UW QSE web site at http : // www .mrf m . org , specifically at http 7] 
//courses . Washington . edu/goodall/MRFM/historical_background . htinl#kay 



REFERENCES 102 

[152] Percival I 1998 Quantum State Diffusion Cambridge University Press. 

[153] Perelomov A 1986 Generalized Coherent States and Their Applications Springer- 
Verlag. 

[154] Perelomov A M 1972 Communications in Mathematical Physics 26, 222-236. 

[155] Peres A and Terno D R 2004 Reviews of Modern Physics 76(1), 93-123. 

[156] Perez-Garcia D, Verstraete F, Wolf M M and Cirac J I 2006 'Matrix product state 
representations'. |arXiv. org:quant-ph/0608197j 

[157] Perez-Garcia D, Verstraete F, Wolf M M and Cirac J I 2007 Quantum Inf. Comput. 
7(5-6), 401-430. 

[158] Powell H F, Eijkelenborg M A V, Irvinel W, Segall D M and Thompson R C 2002 
Journal of Physics B: Atomic, Molecular, and Optical Physics pp. 205-226. 

[159] Preskill J 2006 Course information for Physics 219/Computer Science 219. On-line 
notes: |http : / / www . theory . caltech . edu/people/preskill/ph229j 

[160] Press W H, Flannery B P, Teukolsky S A and Vetterling W T 1994 Numerical Recipes: 
The Art of Scientific Computing, Second Edition, Cambridge University Press Cam- 
bridge. 

[161] Raab A 2000 Chemical Physics Letters 319, 674-8. 
[162] Radchffe J M 1971 Journal of Physics A 4, 313-323. 

[163] Raman S, Qian B, Baker D and Walker R C 2008 IBM Journal of Research and 
Development 52{l/2), 7-18. 

[164] Rewiehski M J 2003 A Trajectory Piecewise-Linear Approach to Model Order Reduc- 
tion of Nonlinear Dynamical Systems PhD thesis MIT. 

[165] Rewiehski M and White J 2006 Linear Algebra and its Applications 415(2-3), 426- 
454. 

[166] Rigo M and Gisin N 1996 Quantum and Semiclassical Optics 8, 255-268. 
[167] Ritov Y 2007 Ann. Statist. 35(6), 2370-2372. 

[168] Rogers L C G and Williams D 2000 Diffusions, Markov Processes and Martingales, 
Volume 2: ltd Calculus Cambridge Mathematical Library. 

[169] Rose M E 1957 Elementary Theory of Angular Momentum John Wiley & Sons. 

[170] Rugar D, Budakian R, Mamin H J and Chui B W 2004 Nature 430, 329-332. 

[171] Rugar D, Mamin H J and Guethner P 1989 Applied Physics Letters 55(25), 2588-90. 

[172] Schilling T A 1996 Geometry of Quantum Mechanics. PhD thesis Pennsylvania State 
University. 

[173] SchoUwock U 2005 Reviews of Modern Physics 77(1), 259-315. 



REFERENCES 103 



[174 
[175; 

[176 

[177; 
[178 
[179 

[180 

[181 

[182; 

[183 
[184 
[185; 
[186 

[187; 

[188 
[189 

[190 



[191 
[192 



Schuch N, Wolf M M, Vollbrecht K G H and Cirac J I 2008 On entropy growth and 
the hardness of simulating time evolution, http://arxiv.org/abs/0801.2078, 

Schurrmacher A 1997 in G Meggle, ed., 'Analyomen 2, Band I: Logic, Epistemology, 
Philosophy of Science' Gesellschaft fiir Analytische Philosophie Walter de Gruyter 
(pub.) Berlin pp. 453-461. 

Shimoda K, Takahasi H and Townes C 1957/06/ Journal of the Physical Society of 
Japan 12(686-700). 

Sidles J A 1991 Applied Physics Letters 58(24), 2854-6. 
Sidles J A 1992 Physical Review Letters 68(8), 1124-7. 

Sidles J A, Garbini J L, Bruland K J, Rugar D, Ziiger O, Hoen S and Yannoni C S 
1995 Reviews of Modern Physics 67(1), 249-265. 

Sidles J A, Garbini J L, Dougherty W M and Chao S H 2003 Proceedings of the IEEE 
91(5), 799-816. 

Sidles J A, Garbini J L and Drobny G P 1992 Review of Scientific Instruments 
63(8), 3881-99. 

Sidles J A and Rugar D 1993 Physical Review Letters 70(22), 3506-9. 
Silverstein J W 1985 The Annals of Probability 13(4), 1364-1368. 
Slater J C 1929 Phys. Rev. 34(10), 1293-1322. 

Spielman D A and Teng S H 2004 J. ACM 51(3), 385-463 (electronic). 

Stinespring W F 1955 Proceedings of the American Mathematical Society 6(2), 211- 
216. 

Tao T 2007 Why global regularity for Navier-Stokes is hard. To appear in a collection. 
What's new - 2007: Open Questions, Expository Articles, and Lecture series from a 
Mathematical Blog. 

Tenenbaum J B, Silva V d and Langford J C 2000 Science 290(5500), 2319-2323. 

Udri§te C 1994 Convex functions and optimization methods on Riemannian mani- 
folds Vol. 297 of Mathematics and its Applications Kluwer Academic Publishers Group 
Dordrecht. 

van den Berg E, Friedlander M P, Hennenfent G, Herrmann F J, Saab R and Ozgur 
Yilmaz 2007 Sparco: A testing framework for sparse reconstruction Technical Report 



TR-2007-2 University of British Columbia, available at http : // www . cs . ubc . ca/labs/ 



scl/sparco/[ 



Vidal G 2004 Physical Review Letters 93(4), 040502 - 1. 

Viollet-Le-Duc E E 1856 Dictionnaire Raisonne de L' Architecture Frangaise du XF 
au XVF Siecle A. Morel, Paris. 



REFERENCES 104 

[193] Vranas P, Blumrich M A, Chen D, Gara A, Giampapa M E, Heidelberger P, Salapura 
V, Sexton J C, Soltz R and Bhanot G 2008 IBM Journal of Research and Development 
52(1/2), 189-198. 

[194] Wakin M and Baraniuk R 2006 Acoustics, Speech and Signal Processing, 2006. 
ICASSP 2006 Proceedings. 2006 IEEE International Conference on 5, V-V. 

[195] Weinberg S 1972 Gravitation and Cosmology: Principles and Applications of the Gen- 
eral Theory of Relativity John Wiley &; Sons. 

[196] Wheeler J and Feynman R 1945 Reviews of Modern Physics 17, 157 - 181. 

[197] Wheeler J and Feynman R 1949 Reviews of Modern Physics 21, 425 - 434. 

[198] Wigner E P 1959 Group Theory and Its Application to the Quantum Mechanics of 
Atomic Spectra Academic Press. 

[199] Wootters W 1998 Phys. Rev. Lett. (USA) 80(10), 2245-8. 

[200] Yau S T 2006 'Perspectives on geometric analysis'. ArXiv preprint math/0602363, to 
appear in Survey in Differential Geometry, Vol. X. 

[201] Zurek W H 2003 Reviews of Modem Physics 75(3), 715-75. 



